TABLE OF CONTENTS
- ABINIT/m_hu
- ABINIT/udens_inglis_hu
- ABINIT/udens_slatercondon_hu
- ABINIT/vee_ylm2jmj_hu
- m_hu/copy_hu
- m_hu/destroy_hu
- m_hu/hu_type
- m_hu/init_hu
- m_hu/print_hu
- m_hu/printvee_hu
- m_hu/reddd
- m_hu/rotatevee_hu
- m_hu/vee2udens_hu
- m_hu/vee2udensatom_hu
- m_hu/vee_ndim2tndim_hu
- m_hu/vee_ndim2tndim_hu_r
- m_hu/vee_slm2ylm_hu
ABINIT/m_hu [ Modules ]
NAME
m_hu
FUNCTION
COPYRIGHT
Copyright (C) 2006-2022 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 .
INPUTS
OUTPUT
SOURCE
19 #if defined HAVE_CONFIG_H 20 #include "config.h" 21 #endif 22 23 24 #include "abi_common.h" 25 26 MODULE m_hu 27 28 use m_abicore 29 30 use defs_basis 31 use m_pawtab, only : pawtab_type 32 33 implicit none 34 35 private 36 37 public :: init_hu 38 public :: copy_hu 39 public :: destroy_hu 40 ! public :: qmc_hu 41 public :: print_hu 42 public :: vee2udens_hu 43 public :: rotatevee_hu 44 public :: printvee_hu 45 public :: vee2udensatom_hu 46 public :: vee_slm2ylm_hu 47 public :: vee_ndim2tndim_hu 48 public :: vee_ndim2tndim_hu_r 49 public :: udens_slatercondon_hu 50 public :: udens_inglis_hu
ABINIT/udens_inglis_hu [ Functions ]
NAME
udens_inglis_hu
FUNCTION
For a given angular momentum l and Slater integrals, give the density density interactions U(m,m') and J(m,m') from Inglis tables in JMJ Basis
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (BA) 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
lcor= angular momentum fk(lcor+1)= Slater integrals
SIDE EFFECTS
NOTES
SOURCE
2250 subroutine udens_inglis_hu(fk,lcor) 2251 2252 use defs_basis 2253 use m_errors 2254 use m_abicore 2255 implicit none 2256 2257 !Arguments --------------------------------------------- 2258 !scalars 2259 integer,intent(in) :: lcor 2260 real(dp), intent(in) :: fk(0:lcor) 2261 2262 !Local variables --------------------------------------- 2263 !scalars 2264 character(len=500) :: message 2265 integer :: m1,m2,kk,tndim 2266 !arrays 2267 real(dp), allocatable :: udens(:,:) 2268 real(dp), allocatable :: jdens(:,:) 2269 real(dp),allocatable :: app(:,:,:) 2270 real(dp),allocatable :: bpp(:,:,:) 2271 real(dp),allocatable :: a2pp(:,:) 2272 real(dp),allocatable :: b2pp(:,:) 2273 2274 !********************************************************************* 2275 tndim=2*(2*lcor+1) 2276 ABI_MALLOC(app,(0:lcor,tndim,tndim)) 2277 ABI_MALLOC(bpp,(0:lcor,tndim,tndim)) 2278 ABI_MALLOC(a2pp,(tndim,tndim)) 2279 ABI_MALLOC(b2pp,(tndim,tndim)) 2280 2281 ABI_MALLOC(udens,(tndim,tndim)) 2282 ABI_MALLOC(jdens,(tndim,tndim)) 2283 udens=zero 2284 jdens=zero 2285 a2pp=zero 2286 b2pp=zero 2287 app=zero 2288 bpp=zero 2289 if(lcor==1) then 2290 app(0,:,:)=one 2291 a2pp(:,:)=RESHAPE((/0._dp,0._dp, 0._dp, 0._dp, 0._dp, 0._dp,& 2292 & 0._dp,0._dp, 0._dp, 0._dp, 0._dp, 0._dp,& 2293 & 0._dp,0._dp, 1._dp,-1._dp,-1._dp, 1._dp,& 2294 & 0._dp,0._dp,-1._dp, 1._dp, 1._dp,-1._dp,& 2295 & 0._dp,0._dp,-1._dp, 1._dp, 1._dp,-1._dp,& 2296 & 0._dp,0._dp, 1._dp,-1._dp,-1._dp, 1._dp/),(/6,6/)) 2297 app(1,:,:)=a2pp(:,:) 2298 app(1,:,:)=app(1,:,:)/25._dp 2299 2300 b2pp(:,:)=RESHAPE((/1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,& 2301 & 0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,& 2302 & 0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,& 2303 & 0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,& 2304 & 0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,& 2305 & 0._dp,0._dp,0._dp,0._dp,0._dp,1._dp /),(/6,6/)) 2306 bpp(0,:,:)=b2pp(:,:) 2307 b2pp(:,:)=RESHAPE((/0._dp,0._dp,1._dp,2._dp,3._dp,4._dp,& 2308 & 0._dp,0._dp,4._dp,3._dp,2._dp,1._dp,& 2309 & 1._dp,4._dp,1._dp,2._dp,2._dp,0._dp,& 2310 & 2._dp,3._dp,2._dp,1._dp,0._dp,2._dp,& 2311 & 3._dp,2._dp,2._dp,0._dp,1._dp,2._dp,& 2312 & 4._dp,1._dp,0._dp,2._dp,2._dp,1._dp /),(/6,6/)) 2313 bpp(1,:,:)=b2pp(:,:) 2314 bpp(1,:,:)=bpp(1,:,:)/25._dp 2315 endif 2316 if(lcor==2) then 2317 app(0,:,:)=one 2318 a2pp(:,:)=RESHAPE((/ 49._dp,-49._dp,-49._dp, 49._dp, 70._dp,-14._dp,-56._dp,-56._dp,-14._dp, 70._dp,& 2319 & -49._dp, 49._dp, 49._dp,-49._dp,-70._dp, 14._dp, 56._dp, 56._dp, 14._dp,-70._dp,& 2320 & -49._dp, 49._dp, 49._dp,-49._dp,-70._dp, 14._dp, 56._dp, 56._dp, 14._dp,-70._dp,& 2321 & 49._dp,-49._dp,-49._dp, 49._dp, 70._dp,-14._dp,-56._dp,-56._dp,-14._dp, 70._dp,& 2322 & 70._dp,-70._dp,-70._dp, 70._dp,100._dp,-20._dp,-80._dp,-80._dp,-20._dp,100._dp,& 2323 & -14._dp, 14._dp, 14._dp,-14._dp,-20._dp, 4._dp, 16._dp, 16._dp, 4._dp,-20._dp,& 2324 & -56._dp, 56._dp, 56._dp,-56._dp,-80._dp, 16._dp, 64._dp, 64._dp, 16._dp,-80._dp,& 2325 & -56._dp, 56._dp, 56._dp,-56._dp,-80._dp, 16._dp, 64._dp, 64._dp, 16._dp,-80._dp,& 2326 & -14._dp, 14._dp, 14._dp,-14._dp,-20._dp, 4._dp, 16._dp, 16._dp, 4._dp,-20._dp,& 2327 & 70._dp,-70._dp,-70._dp, 70._dp,100._dp,-20._dp,-80._dp,-80._dp,-20._dp,100._dp/),(/10,10/)) 2328 app(1,:,:)=a2pp(:,:)/1225._dp 2329 app(2,:,:)=zero 2330 2331 b2pp(:,:)=RESHAPE((/1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,& 2332 & 0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,& 2333 & 0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,& 2334 & 0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,& 2335 & 0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,& 2336 & 0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,& 2337 & 0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,& 2338 & 0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,& 2339 & 0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,& 2340 & 0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp/),(/10,10/)) 2341 bpp(0,:,:)=b2pp(:,:) 2342 b2pp(:,:)=RESHAPE((/ 49._dp, 98._dp, 98._dp, 0._dp, 30._dp, 36._dp, 27._dp, 12._dp, 0._dp, 0._dp,& 2343 & 98._dp, 49._dp, 0._dp, 98._dp, 40._dp, 2._dp, 6._dp, 25._dp, 32._dp, 0._dp,& 2344 & 98._dp, 0._dp, 49._dp, 98._dp, 0._dp, 32._dp, 25._dp, 6._dp, 2._dp, 40._dp,& 2345 & 0._dp, 98._dp, 98._dp, 49._dp, 0._dp, 0._dp, 12._dp, 27._dp, 36._dp, 30._dp,& 2346 & 30._dp, 40._dp, 0._dp, 0._dp,100._dp,120._dp, 60._dp, 0._dp, 0._dp, 0._dp,& 2347 & 36._dp, 2._dp, 32._dp, 0._dp,120._dp, 4._dp, 48._dp,108._dp, 0._dp, 0._dp,& 2348 & 27._dp, 6._dp, 25._dp, 12._dp, 60._dp, 48._dp, 64._dp, 0._dp,108._dp, 0._dp,& 2349 & 12._dp, 25._dp, 6._dp, 27._dp, 0._dp,108._dp, 0._dp, 64._dp, 48._dp, 60._dp,& 2350 & 0._dp, 32._dp, 2._dp, 36._dp, 0._dp, 0._dp,108._dp, 48._dp, 4._dp,120._dp,& 2351 & 0._dp, 0._dp, 40._dp, 30._dp, 0._dp, 0._dp, 0._dp, 60._dp,120._dp,100._dp/),(/10,10/)) 2352 bpp(1,:,:)=b2pp(:,:)/1225._dp 2353 write(std_out,*) "warning: this test is only valid if f4of2_sla=0" 2354 endif 2355 if(lcor==3) then 2356 ! app(0,:,:)=one 2357 ! a2pp(:,:)=RESHAPE((/ 100.0/1225.0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2358 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2359 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2360 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2361 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2362 !& 0 ,0, 0, 0, 0, 100.0/1225.0, 0, 0, 0, 0, 0, 0, 0, 0,& 2363 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2364 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2365 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2366 !& 0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,& 2367 !& 0 ,0, 1,-1,-1, 1, 0, 0, 0, 0, 0, 0, 0, 0,& 2368 !& 0 ,0,-1, 1, 1,-1, 0, 0, 0, 0, 0, 0, 0, 0,& 2369 !& 0 ,0,-1, 1, 1,-1, 0, 0, 0, 0, 0, 0, 0, 0,& 2370 !& 0 ,0, 1,-1,-1, 1, 0, 0, 0, 0, 0, 0, 0, 0,&/),(/14,14/)) 2371 b2pp(:,:)=RESHAPE((/ 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2372 & 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2373 & 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2374 & 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2375 & 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2376 & 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2377 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2378 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2379 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2380 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0,& 2381 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0,& 2382 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0,& 2383 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0,& 2384 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0/),(/14,14/)) 2385 bpp(0,:,:)=b2pp(:,:) 2386 b2pp(:,:)=RESHAPE((/100.0, 120.0, 60.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2387 & 120.0, 4.0, 48.0, 108.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2388 & 60.0, 48.0, 64.0, 0.0, 108.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2389 & 0.0, 108.0, 0.0, 64.0, 48.0, 60.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2390 & 0.0, 0.0, 108.0, 48.0, 4.0, 120.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2391 & 0.0, 0.0, 0.0, 60.0, 120.0, 100.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2392 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2393 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2394 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2395 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2396 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2397 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2398 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,& 2399 & 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0/),(/14,14/)) 2400 bpp(1,:,:)=b2pp(:,:)/1225._dp 2401 write(std_out,*) "warning: only 5/2 5/2 elements are given and only for exchange and F2 !" 2402 ! app(1,:,:)=a2pp(:,:) 2403 ! app(1,:,:)=app(1,:,:)/25 2404 ! 2405 ! bpp(0,:,:)=one 2406 ! b2pp(:,:)=RESHAPE((/0,0,1,2,3,4,& 2407 !& 0,0,4,3,2,1,& 2408 !& 1,4,1,2,2,0,& 2409 !& 2,3,2,1,0,2,& 2410 !& 3,2,2,0,1,2,& 2411 !& 4,1,0,2,2,1 /),(/6,6/)) 2412 ! bpp(1,:,:)=b2pp(:,:) 2413 ! bpp(1,:,:)=bpp(1,:,:)/25 2414 endif 2415 2416 do kk=0,lcor 2417 do m1=1,tndim 2418 do m2=1,tndim 2419 !write(6,*) kk,m1,m2 2420 !write(6,*) "--",fk(kk),aklmlmp(kk,m1,m2) 2421 !write(6,*) "--",fk(kk),bklmlmp(kk,m1,m2) 2422 udens(m1,m2)=udens(m1,m2)+fk(kk)*app(kk,m1,m2) 2423 jdens(m1,m2)=jdens(m1,m2)+fk(kk)*bpp(kk,m1,m2) 2424 enddo 2425 enddo 2426 enddo 2427 write(message,'(2x,a,3x,14f10.4)') " Direct Interaction Matrix from Inglis tables (in the JMJ basis) " 2428 call wrtout(std_out, message,'COLL') 2429 write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,tndim,1) 2430 call wrtout(std_out, message,'COLL') 2431 do m1=1,tndim 2432 write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2),m2=1,tndim,1) 2433 call wrtout(std_out, message,'COLL') 2434 enddo 2435 2436 write(message,'(a,2x,a,3x,14f10.4)') ch10," Exchange Interaction Matrix from Inglis tables (in the JMJ basis) " 2437 call wrtout(std_out, message,'COLL') 2438 write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,tndim,1) 2439 call wrtout(std_out, message,'COLL') 2440 do m1=1,tndim 2441 write(message,'(2x,i4,3x,14f10.4)') m1,(jdens(m1,m2),m2=1,tndim,1) 2442 call wrtout(std_out, message,'COLL') 2443 enddo 2444 2445 write(message,'(a,2x,a,3x,14f10.4)') ch10, " Density Density interactions from Inglis tables (in the JMJ basis) " 2446 call wrtout(std_out, message,'COLL') 2447 write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,tndim,1) 2448 call wrtout(std_out, message,'COLL') 2449 do m1=1,tndim 2450 write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2)-jdens(m1,m2),m2=1,tndim,1) 2451 call wrtout(std_out, message,'COLL') 2452 enddo 2453 2454 2455 ABI_FREE(jdens) 2456 ABI_FREE(udens) 2457 ABI_FREE(app) 2458 ABI_FREE(bpp) 2459 ABI_FREE(a2pp) 2460 ABI_FREE(b2pp) 2461 2462 end subroutine udens_inglis_hu
ABINIT/udens_slatercondon_hu [ Functions ]
NAME
udens_slatercondon_hu
FUNCTION
For a given angular momentum l and Slater integrals, give the density density interactions U(m,m') and J(m,m') from Slater and Condon tables
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (BA) 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
lcor= angular momentum fk(lcor+1)= Slater integrals
SIDE EFFECTS
NOTES
SOURCE
1981 subroutine udens_slatercondon_hu(fk,lcor) 1982 1983 use defs_basis 1984 use m_errors 1985 use m_abicore 1986 implicit none 1987 1988 !Arguments --------------------------------------------- 1989 !scalars 1990 integer,intent(in) :: lcor 1991 real(dp), intent(in) :: fk(0:lcor) 1992 1993 !Local variables --------------------------------------- 1994 !scalars 1995 character(len=500) :: message 1996 integer :: m1,m2,kk 1997 !arrays 1998 real(dp), allocatable :: aklmlmp(:,:,:) 1999 real(dp), allocatable :: bklmlmp(:,:,:) 2000 real(dp), allocatable :: udens(:,:) 2001 real(dp), allocatable :: jdens(:,:) 2002 2003 !********************************************************************* 2004 ABI_MALLOC(aklmlmp,(0:lcor,-lcor:lcor,-lcor:lcor)) ! k,m,m' 2005 ABI_MALLOC(bklmlmp,(0:lcor,-lcor:lcor,-lcor:lcor)) ! k,m,m' 2006 ABI_MALLOC(udens,(-lcor:lcor,-lcor:lcor)) ! m,m' 2007 ABI_MALLOC(jdens,(-lcor:lcor,-lcor:lcor)) ! m,m' 2008 ! k=2*(lcor) 2009 aklmlmp=zero 2010 bklmlmp=zero 2011 udens=zero 2012 jdens=zero 2013 if(lcor==0) then 2014 aklmlmp(0,0,0)=1 2015 ! 2016 bklmlmp(0,0,0)=0 2017 else if(lcor==1) then 2018 aklmlmp(0, :, :)=1 2019 aklmlmp(1, 1, 1)= one/25._dp 2020 aklmlmp(1,-1,-1)= one/25._dp 2021 aklmlmp(1,-1, 1)= one/25._dp 2022 aklmlmp(1, 1,-1)= one/25._dp 2023 aklmlmp(1, 1, 0)=-two/25._dp 2024 aklmlmp(1,-1, 0)=-two/25._dp 2025 aklmlmp(1, 0,-1)=-two/25._dp 2026 aklmlmp(1, 0, 1)=-two/25._dp 2027 aklmlmp(1, 0, 0)= four/25._dp 2028 ! 2029 bklmlmp(0, 1, 1)= one 2030 bklmlmp(0,-1,-1)= one 2031 bklmlmp(0, 0, 0)= one 2032 bklmlmp(1, 1, 1)= one/25._dp 2033 bklmlmp(1,-1,-1)= one/25._dp 2034 bklmlmp(1,-1,+1)= six/25._dp 2035 bklmlmp(1,+1,-1)= six/25._dp 2036 bklmlmp(1,+1, 0)= three/25._dp 2037 bklmlmp(1,-1, 0)= three/25._dp 2038 bklmlmp(1, 0,-1)= three/25._dp 2039 bklmlmp(1, 0, 1)= three/25._dp 2040 bklmlmp(1, 0, 0)= four/25._dp 2041 else if(lcor==2) then 2042 aklmlmp(0, :, :)=1 2043 aklmlmp(1, 2, 2)= four / 49._dp 2044 aklmlmp(1, 2, 1)= -two / 49._dp 2045 aklmlmp(1, 2, 0)= -four / 49._dp 2046 aklmlmp(1, 2,-1)= -two / 49._dp 2047 aklmlmp(1, 2,-2)= four / 49._dp 2048 aklmlmp(1, 1, 2)= -two / 49._dp 2049 aklmlmp(1, 1, 1)= one / 49._dp 2050 aklmlmp(1, 1, 0)= two / 49._dp 2051 aklmlmp(1, 1,-1)= one / 49._dp 2052 aklmlmp(1, 1,-2)= -two / 49._dp 2053 aklmlmp(1, 0, 2)= -four / 49._dp 2054 aklmlmp(1, 0, 1)= two / 49._dp 2055 aklmlmp(1, 0, 0)= four / 49._dp 2056 aklmlmp(1, 0,-1)= two / 49._dp 2057 aklmlmp(1, 0,-2)= -four / 49._dp 2058 aklmlmp(1,-1, 2)= -two / 49._dp 2059 aklmlmp(1,-1, 1)= one / 49._dp 2060 aklmlmp(1,-1, 0)= two / 49._dp 2061 aklmlmp(1,-1,-1)= one / 49._dp 2062 aklmlmp(1,-1,-2)= -two / 49._dp 2063 aklmlmp(1,-2, 2)= four / 49._dp 2064 aklmlmp(1,-2, 1)= -two / 49._dp 2065 aklmlmp(1,-2, 0)= -four / 49._dp 2066 aklmlmp(1,-2,-1)= -two / 49._dp 2067 aklmlmp(1,-2,-2)= four / 49._dp 2068 2069 aklmlmp(2, 2, 2)= one / 441._dp 2070 aklmlmp(2, 2, 1)= -four / 441._dp 2071 aklmlmp(2, 2, 0)= six / 441._dp 2072 aklmlmp(2, 2,-1)= -four / 441._dp 2073 aklmlmp(2, 2,-2)= one / 441._dp 2074 aklmlmp(2, 1, 2)= -four / 441._dp 2075 aklmlmp(2, 1, 1)= 16._dp / 441._dp 2076 aklmlmp(2, 1, 0)= -24._dp / 441._dp 2077 aklmlmp(2, 1,-1)= 16._dp / 441._dp 2078 aklmlmp(2, 1,-2)= -four / 441._dp 2079 aklmlmp(2, 0, 2)= six / 441._dp 2080 aklmlmp(2, 0, 1)= -24._dp / 441._dp 2081 aklmlmp(2, 0, 0)= 36._dp / 441._dp 2082 aklmlmp(2, 0,-1)= -24._dp / 441._dp 2083 aklmlmp(2, 0,-2)= six / 441._dp 2084 aklmlmp(2,-1, 2)= -four / 441._dp 2085 aklmlmp(2,-1, 1)= 16._dp / 441._dp 2086 aklmlmp(2,-1, 0)= -24._dp / 441._dp 2087 aklmlmp(2,-1,-1)= 16._dp / 441._dp 2088 aklmlmp(2,-1,-2)= -four / 441._dp 2089 aklmlmp(2,-2, 2)= one / 441._dp 2090 aklmlmp(2,-2, 1)= -four / 441._dp 2091 aklmlmp(2,-2, 0)= six / 441._dp 2092 aklmlmp(2,-2,-1)= -four / 441._dp 2093 aklmlmp(2,-2,-2)= one / 441._dp 2094 !do m1=lcor,-lcor,-1 2095 ! do m2=lcor,-lcor,-1 2096 ! write(6,*) m1,m2,aklmlmp(2,m1,m2) 2097 ! enddo 2098 !enddo 2099 2100 !write(message,'(2x,a,3x,14f10.4)') " Slater aklmlmp(2,:,:)" 2101 !call wrtout(std_out, message,'COLL') 2102 !write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1) 2103 !call wrtout(std_out, message,'COLL') 2104 !do m1=-lcor,lcor,1 2105 ! write(message,'(2x,i4,3x,14f10.4)') m1,(aklmlmp(2,m1,m2),m2=-lcor,lcor,1) 2106 ! call wrtout(std_out, message,'COLL') 2107 !enddo 2108 2109 bklmlmp(0, 2, 2)=1 2110 bklmlmp(0, 1, 1)=1 2111 bklmlmp(0, 0, 0)=1 2112 bklmlmp(0,-1,-1)=1 2113 bklmlmp(0,-2,-2)=1 2114 bklmlmp(1, 2, 2)= four / 49._dp 2115 bklmlmp(1, 2, 1)= six / 49._dp 2116 bklmlmp(1, 2, 0)= four / 49._dp 2117 bklmlmp(1, 2,-1)= zero 2118 bklmlmp(1, 2,-2)= zero 2119 bklmlmp(1, 1, 2)= six / 49._dp 2120 bklmlmp(1, 1, 1)= one / 49._dp 2121 bklmlmp(1, 1, 0)= one / 49._dp 2122 bklmlmp(1, 1,-1)= six / 49._dp 2123 bklmlmp(1, 1,-2)= zero 2124 bklmlmp(1, 0, 2)= four / 49._dp 2125 bklmlmp(1, 0, 1)= one / 49._dp 2126 bklmlmp(1, 0, 0)= four / 49._dp 2127 bklmlmp(1, 0,-1)= one / 49._dp 2128 bklmlmp(1, 0,-2)= four / 49._dp 2129 bklmlmp(1,-1, 2)= zero 2130 bklmlmp(1,-1, 1)= six / 49._dp 2131 bklmlmp(1,-1, 0)= one / 49._dp 2132 bklmlmp(1,-1,-1)= one / 49._dp 2133 bklmlmp(1,-1,-2)= six / 49._dp 2134 bklmlmp(1,-2, 2)= zero 2135 bklmlmp(1,-2, 1)= zero 2136 bklmlmp(1,-2, 0)= four / 49._dp 2137 bklmlmp(1,-2,-1)= six / 49._dp 2138 bklmlmp(1,-2,-2)= four / 49._dp 2139 2140 bklmlmp(2, 2, 2)= one / 441._dp 2141 bklmlmp(2, 2, 1)= five / 441._dp 2142 bklmlmp(2, 2, 0)= 15._dp / 441._dp 2143 bklmlmp(2, 2,-1)= 35._dp / 441._dp 2144 bklmlmp(2, 2,-2)= 70._dp / 441._dp 2145 bklmlmp(2, 1, 2)= five / 441._dp 2146 bklmlmp(2, 1, 1)= 16._dp / 441._dp 2147 bklmlmp(2, 1, 0)= 30._dp / 441._dp 2148 bklmlmp(2, 1,-1)= 40._dp / 441._dp 2149 bklmlmp(2, 1,-2)= 35._dp / 441._dp 2150 bklmlmp(2, 0, 2)= 15._dp / 441._dp 2151 bklmlmp(2, 0, 1)= 30._dp / 441._dp 2152 bklmlmp(2, 0, 0)= 36._dp / 441._dp 2153 bklmlmp(2, 0,-1)= 30._dp / 441._dp 2154 bklmlmp(2, 0,-2)= 15._dp / 441._dp 2155 bklmlmp(2,-1, 2)= 35._dp / 441._dp 2156 bklmlmp(2,-1, 1)= 40._dp / 441._dp 2157 bklmlmp(2,-1, 0)= 30._dp / 441._dp 2158 bklmlmp(2,-1,-1)= 16._dp / 441._dp 2159 bklmlmp(2,-1,-2)= five / 441._dp 2160 bklmlmp(2,-2, 2)= 70._dp / 441._dp 2161 bklmlmp(2,-2, 1)= 35._dp / 441._dp 2162 bklmlmp(2,-2, 0)= 15._dp / 441._dp 2163 bklmlmp(2,-2,-1)= five / 441._dp 2164 bklmlmp(2,-2,-2)= one / 441._dp 2165 2166 !write(message,'(2x,a,3x,14f10.4)') " Slater bklmlmp(2,:,:)" 2167 !call wrtout(std_out, message,'COLL') 2168 !write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1) 2169 !call wrtout(std_out, message,'COLL') 2170 !do m1=-lcor,lcor,1 2171 ! write(message,'(2x,i4,3x,14f10.4)') m1,(bklmlmp(2,m1,m2),m2=-lcor,lcor,1) 2172 ! call wrtout(std_out, message,'COLL') 2173 !enddo 2174 ! if f4of2_sla: these data agree with the explicit calculation 2175 else if(lcor==3) then 2176 endif 2177 2178 do kk=0,lcor 2179 do m1=-lcor,lcor,1 2180 do m2=-lcor,lcor,1 2181 !write(6,*) kk,m1,m2 2182 !write(6,*) "--",fk(kk),aklmlmp(kk,m1,m2) 2183 !write(6,*) "--",fk(kk),bklmlmp(kk,m1,m2) 2184 udens(m1,m2)=udens(m1,m2)+fk(kk)*aklmlmp(kk,m1,m2) 2185 jdens(m1,m2)=jdens(m1,m2)+fk(kk)*bklmlmp(kk,m1,m2) 2186 enddo 2187 enddo 2188 enddo 2189 write(message,'(2x,a,3x,14f10.4)') " Direct Interaction Matrix from Slater tables (in the Ylm basis) " 2190 call wrtout(std_out, message,'COLL') 2191 write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1) 2192 call wrtout(std_out, message,'COLL') 2193 do m1=-lcor,lcor,1 2194 write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2),m2=-lcor,lcor,1) 2195 call wrtout(std_out, message,'COLL') 2196 enddo 2197 2198 write(message,'(a,2x,a,3x,14f10.4)') ch10," Exchange Interaction Matrix from Slater tables (in the Ylm basis) " 2199 call wrtout(std_out, message,'COLL') 2200 write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1) 2201 call wrtout(std_out, message,'COLL') 2202 do m1=-lcor,lcor,1 2203 write(message,'(2x,i4,3x,14f10.4)') m1,(jdens(m1,m2),m2=-lcor,lcor,1) 2204 call wrtout(std_out, message,'COLL') 2205 enddo 2206 2207 write(message,'(a,2x,a,3x,14f10.4)') ch10," Density density Interaction Matrix from Slater tables (in the Ylm basis) " 2208 call wrtout(std_out, message,'COLL') 2209 write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1) 2210 call wrtout(std_out, message,'COLL') 2211 do m1=-lcor,lcor,1 2212 write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2)-jdens(m1,m2),m2=-lcor,lcor,1) 2213 call wrtout(std_out, message,'COLL') 2214 enddo 2215 2216 ABI_FREE(jdens) 2217 ABI_FREE(udens) 2218 ABI_FREE(aklmlmp) 2219 ABI_FREE(bklmlmp) 2220 2221 end subroutine udens_slatercondon_hu
ABINIT/vee_ylm2jmj_hu [ Functions ]
NAME
vee_ylm2jmj_hu
FUNCTION
For a given angular momentum lcor, change a matrix of dimension [2(2*lcor+1)]**4 from the Ylm basis to the J,M_J basis if option==1
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (BA) 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
lcor= angular momentum option= 1 matrix in |l,s,m_l,m_s> basis is changed into |l,s,j,m_j> basis 2 matrix in |l,s,j,m_j> basis is changed into |l,s,m_l,m_s> basis
SIDE EFFECTS
mat_mlms= Input/Ouput matrix in the Ylm basis, size of the matrix is (2*lcor+1,2*lcor+1,ndij) mat_jmj= Input/Output matrix in the J,M_J basis, size is 2*(2*lcor+1),2*(2*lcor+1)
NOTES
usefull only in ndij==4
SOURCE
1828 subroutine vee_ylm2jmj_hu(lcor,mat_inp_c,mat_out_c,option) 1829 1830 use defs_basis 1831 use m_errors 1832 use m_abicore 1833 implicit none 1834 1835 !Arguments --------------------------------------------- 1836 !scalars 1837 integer,intent(in) :: lcor,option 1838 !arrays 1839 complex(dpc),intent(in) :: mat_inp_c(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1)) 1840 complex(dpc),intent(out) :: mat_out_c(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1)) 1841 1842 !Local variables --------------------------------------- 1843 !scalars 1844 integer :: ii,im,jc1,jj,jm,ll,ml1,ms1,gg,hh,gm,hm 1845 real(dp),parameter :: invsqrt2=one/sqrt2 1846 real(dp) :: invsqrt2lp1,xj,xmj 1847 complex(dpc) :: tmp2 1848 character(len=500) :: message 1849 !arrays 1850 integer, allocatable :: ind_msml(:,:) 1851 complex(dpc),allocatable :: mlms2jmj(:,:) 1852 1853 !********************************************************************* 1854 1855 if (option/=1.and.option/=2) then 1856 message=' option=/1 and =/2 !' 1857 ABI_BUG(message) 1858 end if 1859 1860 if(option==1) then 1861 write(message,'(3a)') ch10,& 1862 & "matrix in |l,s,m_l,m_s> basis is changed into |l,s,j,m_j> basis" 1863 call wrtout(std_out,message) 1864 else if(option==2) then 1865 write(message,'(3a)') ch10,& 1866 & "matrix in |l,s,j,m_j> basis is changed into |l,s,m_l,m_s> basis" 1867 call wrtout(std_out,message) 1868 end if 1869 1870 !--------------- Built indices + allocations 1871 ll=lcor 1872 ABI_MALLOC(mlms2jmj,(2*(2*ll+1),2*(2*ll+1))) 1873 mlms2jmj=czero 1874 ABI_MALLOC(ind_msml,(2,-ll:ll)) 1875 mlms2jmj=czero 1876 jc1=0 1877 do ms1=1,2 1878 do ml1=-ll,ll 1879 jc1=jc1+1 1880 ind_msml(ms1,ml1)=jc1 1881 end do 1882 end do 1883 1884 !--------------- built mlms2jmj 1885 !do jj=ll,ll+1 ! the physical value of j are ll-0.5,ll+0.5 1886 !xj(jj)=jj-0.5 1887 if(ll==0)then 1888 message=' ll should not be equal to zero !' 1889 ABI_BUG(message) 1890 end if 1891 jc1=0 1892 invsqrt2lp1=one/sqrt(float(2*lcor+1)) 1893 do jj=ll,ll+1 1894 xj=float(jj)-half ! xj is in {ll-0.5, ll+0.5} 1895 do jm=-jj,jj-1 1896 xmj=float(jm)+half ! xmj is in {-xj,xj} 1897 jc1=jc1+1 ! Global index for JMJ 1898 if(nint(xj+0.5)==ll+1) then ! if xj=ll+0.5 1899 if(nint(xmj+0.5)==ll+1) then 1900 mlms2jmj(ind_msml(2,ll),jc1)=1.0 ! J=L+0.5 and m_J=L+0.5 1901 else if(nint(xmj-0.5)==-ll-1) then 1902 mlms2jmj(ind_msml(1,-ll),jc1)=1.0 ! J=L+0.5 and m_J=-L-0.5 1903 else 1904 mlms2jmj(ind_msml(2,nint(xmj-0.5)),jc1)=invsqrt2lp1*(sqrt(float(ll)+xmj+0.5)) 1905 mlms2jmj(ind_msml(1,nint(xmj+0.5)),jc1)=invsqrt2lp1*(sqrt(float(ll)-xmj+0.5)) 1906 end if 1907 end if 1908 if(nint(xj+0.5)==ll) then ! if xj=ll-0.5 1909 mlms2jmj(ind_msml(1,nint(xmj+0.5)),jc1)=invsqrt2lp1*(sqrt(float(ll)+xmj+0.5)) 1910 mlms2jmj(ind_msml(2,nint(xmj-0.5)),jc1)=-invsqrt2lp1*(sqrt(float(ll)-xmj+0.5)) 1911 end if 1912 end do 1913 end do 1914 write(message,'(3a)') ch10,"Matrix to go from |M_L,M_S> to |J,M_J>" 1915 call wrtout(std_out,message,"COLL") 1916 do im=1,2*(ll*2+1) 1917 write(message,'(12(1x,18(1x,"(",f7.3,",",f7.3,")")))') (mlms2jmj(im,jm),jm=1,2*(ll*2+1)) 1918 call wrtout(std_out,message,"COLL") 1919 end do 1920 1921 !--------------- compute change of basis 1922 do jm=1,2*(2*ll+1) 1923 do im=1,2*(2*ll+1) 1924 do hm=1,2*(2*ll+1) 1925 do gm=1,2*(2*ll+1) 1926 tmp2=czero 1927 do gg=1,2*(2*ll+1) 1928 do hh=1,2*(2*ll+1) 1929 do ii=1,2*(2*ll+1) 1930 do jj=1,2*(2*ll+1) 1931 if(option==1) then 1932 tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*CONJG(mlms2jmj(ii,im))*(mlms2jmj(jj,jm))& 1933 & *CONJG(mlms2jmj(gg,gm))*(mlms2jmj(hh,hm)) 1934 else if(option==2) then 1935 ! tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*CONJG(mlms2jmj(ii,im))*(mlms2jmj(jj,jm))& ! inv=t* 1936 !& *CONJG(mlms2jmj(gg,gm))*(mlms2jmj(hh,hm)) ! inv=t* 1937 tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*(mlms2jmj(im,ii))*CONJG(mlms2jmj(jm,jj))& ! inv=t* 1938 & *(mlms2jmj(gm,gg))*CONJG(mlms2jmj(hm,hh)) ! inv=t* 1939 end if 1940 end do 1941 end do 1942 end do 1943 end do 1944 mat_out_c(gm,im,hm,jm)=tmp2 1945 end do 1946 end do 1947 end do 1948 end do 1949 ABI_FREE(mlms2jmj) 1950 ABI_FREE(ind_msml) 1951 1952 end subroutine vee_ylm2jmj_hu
m_hu/copy_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
copy_hu
FUNCTION
Allocate variables used in type hu_type.
INPUTS
hu <type(hu_type)>= U interaction OUTPUTS hu_new <type(hu_type)>= U interaction
SOURCE
312 subroutine copy_hu(ntypat,hu,hu_new) 313 314 use defs_basis 315 implicit none 316 317 !Arguments ------------------------------------ 318 !type 319 integer, intent(in) :: ntypat 320 type(hu_type), intent(in) :: hu(ntypat) 321 type(hu_type), intent(out) :: hu_new(ntypat) 322 !Local variables ------------------------------------ 323 integer :: itypat,ndim 324 !************************************************************************ 325 326 do itypat=1,ntypat 327 hu_new(itypat)%lpawu = hu(itypat)%lpawu 328 hu_new(itypat)%jmjbasis = hu(itypat)%jmjbasis 329 hu_new(itypat)%upawu = hu(itypat)%upawu 330 hu_new(itypat)%jpawu = hu(itypat)%jpawu 331 hu_new(itypat)%f2_sla = hu(itypat)%f2_sla 332 hu_new(itypat)%f4of2_sla = hu(itypat)%f4of2_sla 333 hu_new(itypat)%f6of2_sla = hu(itypat)%f6of2_sla 334 hu_new(itypat)%jpawu_zero = hu(itypat)%jpawu_zero 335 ndim=2*hu_new(itypat)%lpawu+1 336 ABI_MALLOC(hu_new(itypat)%uqmc,(ndim*(2*ndim-1))) 337 ABI_MALLOC(hu_new(itypat)%udens,(2*ndim,2*ndim)) 338 ABI_MALLOC(hu_new(itypat)%vee,(ndim,ndim,ndim,ndim)) 339 hu_new(itypat)%vee = hu(itypat)%vee 340 hu_new(itypat)%udens = hu(itypat)%udens 341 hu_new(itypat)%uqmc = hu(itypat)%uqmc 342 enddo ! itypat 343 344 end subroutine copy_hu
m_hu/destroy_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
destroy_hu
FUNCTION
deallocate hu
INPUTS
ntypat = number of species hu <type(hu_type)> = data for the interaction in DMFT.
OUTPUT
SOURCE
362 subroutine destroy_hu(hu,ntypat,t2g,x2my2d) 363 364 use defs_basis 365 use m_crystal, only : crystal_t 366 implicit none 367 368 !Arguments ------------------------------------ 369 !scalars 370 integer, intent(in) :: ntypat 371 type(hu_type),intent(inout) :: hu(ntypat) 372 integer, intent(in) :: t2g,x2my2d 373 374 !Local variables------------------------------- 375 integer :: itypat 376 377 ! ********************************************************************* 378 379 do itypat=1,ntypat 380 ! if ( allocated(hu(itypat)%vee) ) deallocate(hu(itypat)%vee) 381 if ( allocated(hu(itypat)%uqmc) ) then 382 ABI_FREE(hu(itypat)%uqmc) 383 end if 384 if ( allocated(hu(itypat)%udens) ) then 385 ABI_FREE(hu(itypat)%udens) 386 end if 387 if(t2g==1.and.hu(itypat)%lpawu==1) then 388 ABI_FREE(hu(itypat)%vee) 389 else if(x2my2d==1.and.hu(itypat)%lpawu==0) then 390 ABI_FREE(hu(itypat)%vee) 391 else 392 hu(itypat)%vee => null() 393 endif 394 enddo 395 396 end subroutine destroy_hu
m_hu/hu_type [ Types ]
NAME
hu_type
FUNCTION
This structured datatype contains interaction matrices for the correlated subspace
SOURCE
64 type, public :: hu_type ! for each typat 65 66 integer :: lpawu 67 68 logical :: jmjbasis 69 70 real(dp) :: upawu ! => upaw 71 72 real(dp) :: jpawu ! => jpaw 73 74 real(dp) :: f2_sla ! => f2_sla 75 76 real(dp) :: f4of2_sla ! => f4of2_sla 77 78 real(dp) :: f6of2_sla ! => f6of2_sla 79 80 logical :: jpawu_zero ! true if all jpawu are zero 81 ! false if one of the jpaw is not zero 82 83 real(dp), pointer :: vee(:,:,:,:) => null() ! => vee 84 85 real(dp), allocatable :: uqmc(:) 86 87 real(dp), allocatable :: udens(:,:) 88 89 end type hu_type 90 91 !---------------------------------------------------------------------- 92 93 94 CONTAINS !========================================================================================
m_hu/init_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
init_hu
FUNCTION
Allocate variables used in type hu_type.
INPUTS
cryst_struc <type(crystal_t)>=crystal structure data pawtab <type(pawtab)>=paw related data OUTPUTS hu <type(hu_type)>= U interaction
SOURCE
113 subroutine init_hu(cryst_struc,pawtab,hu,t2g,x2my2d) 114 115 use defs_basis 116 use m_crystal, only : crystal_t 117 implicit none 118 119 !Arguments ------------------------------------ 120 !type 121 type(crystal_t),intent(in) :: cryst_struc 122 type(pawtab_type), target, intent(in) :: pawtab(cryst_struc%ntypat) 123 type(hu_type), intent(inout) :: hu(cryst_struc%ntypat) 124 integer, intent(in) :: t2g,x2my2d 125 !Local variables ------------------------------------ 126 integer :: itypat,i,ij,ij1,ij2,j,lpawu,ms,ms1,m,m1,ndim 127 integer :: ns,ns1,n,n1 128 integer, allocatable :: xij(:,:) 129 real(dp) :: xtemp 130 character(len=500) :: message 131 !************************************************************************ 132 write(message,'(2a)') ch10," == Compute Interactions for DMFT" 133 call wrtout(std_out,message,'COLL') 134 135 xtemp=zero 136 137 ! ==================================== 138 ! Compute hu(iatom)%uqmc from vee 139 ! ==================================== 140 hu(1)%jpawu_zero=.true. 141 do itypat=1,cryst_struc%ntypat 142 hu(itypat)%lpawu=pawtab(itypat)%lpawu 143 hu(itypat)%jmjbasis=.false. 144 if(t2g==1.and.hu(itypat)%lpawu==2) hu(itypat)%lpawu=1 145 if(x2my2d==1.and.hu(itypat)%lpawu==2) hu(itypat)%lpawu=0 146 lpawu=hu(itypat)%lpawu 147 if(lpawu.ne.-1) then 148 hu(itypat)%upawu=pawtab(itypat)%upawu 149 hu(itypat)%jpawu=pawtab(itypat)%jpawu 150 151 ! The if below are necessary for an unknown reason with the NAG compiler. 152 if(pawtab(itypat)%f4of2_sla>=0) then 153 hu(itypat)%f4of2_sla = pawtab(itypat)%f4of2_sla 154 else if(pawtab(itypat)%f4of2_sla==0) then 155 hu(itypat)%f4of2_sla = 0_dp 156 else if(pawtab(itypat)%f4of2_sla<0) then 157 hu(itypat)%f4of2_sla = pawtab(itypat)%f4of2_sla 158 endif 159 if(pawtab(itypat)%f6of2_sla>=0) then 160 hu(itypat)%f6of2_sla = pawtab(itypat)%f6of2_sla 161 else if(pawtab(itypat)%f6of2_sla==0) then 162 hu(itypat)%f6of2_sla = 0_dp 163 else if(pawtab(itypat)%f6of2_sla<0) then 164 hu(itypat)%f6of2_sla = pawtab(itypat)%f6of2_sla 165 endif 166 167 ! This is copied from pawpuxinit: it would be better not to duplicate 168 ! these lines. 169 if(lpawu==0) then 170 hu(itypat)%f2_sla=zero 171 else if(lpawu==1) then 172 hu(itypat)%f2_sla=pawtab(itypat)%jpawu*5._dp 173 else if(lpawu==2) then 174 hu(itypat)%f2_sla=pawtab(itypat)%jpawu*14._dp/(One+hu(itypat)%f4of2_sla) 175 else if(lpawu==3) then 176 hu(itypat)%f2_sla=pawtab(itypat)%jpawu*6435._dp/(286._dp+& 177 & 195._dp*hu(itypat)%f4of2_sla+250._dp*hu(itypat)%f6of2_sla) 178 endif 179 180 if(hu(itypat)%jpawu>tol4) hu(1)%jpawu_zero=.false. 181 ndim=2*lpawu+1 182 ! ndim1=2*hu(itypat)%lpawu+1 183 write(message,'(2a,i4)') ch10,' -------> For Correlated Species', itypat 184 call wrtout(std_out, message,'COLL') 185 ! allocate(hu(itypat)%vee(ndim,ndim,ndim,ndim)) 186 ABI_MALLOC(hu(itypat)%uqmc,(ndim*(2*ndim-1))) 187 ABI_MALLOC(hu(itypat)%udens,(2*ndim,2*ndim)) 188 ABI_MALLOC(xij,(2*ndim,2*ndim)) 189 if(t2g==0.and.x2my2d==0) then 190 hu(itypat)%vee => pawtab(itypat)%vee 191 ! t2g case begin 192 else if(t2g==1.and.hu(itypat)%lpawu==1) then 193 ABI_MALLOC(hu(itypat)%vee,(ndim,ndim,ndim,ndim)) 194 n=0 195 do m=1,5 196 if((m/=1.and.m/=2.and.m/=4)) cycle 197 n=n+1 198 ns=0 199 do ms=1,5 200 if((ms/=1.and.ms/=2.and.ms/=4)) cycle 201 ns=ns+1 202 n1=0 203 do m1=1,5 204 if((m1/=1.and.m1/=2.and.m1/=4)) cycle 205 n1=n1+1 206 ns1=0 207 do ms1=1,5 208 if((ms1/=1.and.ms1/=2.and.ms1/=4)) cycle 209 ns1=ns1+1 210 hu(itypat)%vee(n,ns,n1,ns1)=pawtab(itypat)%vee(m,ms,m1,ms1) 211 enddo 212 enddo 213 enddo 214 enddo 215 ! t2g case end 216 ! x2my2d case begin 217 else if(x2my2d==1.and.hu(itypat)%lpawu==0) then 218 ABI_MALLOC(hu(itypat)%vee,(ndim,ndim,ndim,ndim)) 219 hu(itypat)%vee(1,1,1,1)=pawtab(itypat)%upawu 220 endif 221 ! x2my2d case end 222 223 hu(itypat)%udens=zero 224 ij=0 225 do ms=1,2*ndim-1 226 xij(ms,ms)=0 227 do ms1=ms+1,2*ndim 228 ij=ij+1 229 xij(ms,ms1)=ij 230 xij(ms1,ms)=ij 231 if(ms<=ndim.and.ms1>ndim) then 232 m1 = ms1 - ndim 233 m = ms 234 hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1) 235 hu(itypat)%udens(ms,ms1)= hu(itypat)%vee(m,m1,m,m1) 236 hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1) 237 else if(ms<=ndim.and.ms1<=ndim) then 238 m1 = ms1 239 m = ms 240 hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m) 241 hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij) 242 hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1) 243 else 244 m1 = ms1 - ndim 245 m = ms - ndim 246 hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m) 247 hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij) 248 hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1) 249 endif 250 enddo 251 enddo 252 xij(2*ndim,2*ndim)=0 253 write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation " 254 call wrtout(std_out, message,'COLL') 255 write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim) 256 call wrtout(std_out, message,'COLL') 257 ! xtemp1b=0.d0 258 ! ==================================== 259 ! Print hu(iatom)%uqmc 260 ! ==================================== 261 ij1=-10 262 ij2=-10 263 ij=0 264 do i=1,2*ndim 265 do j=i+1,2*ndim 266 ij=ij+1 267 if(j==i+1) ij1=ij 268 if(j==2*ndim) ij2=ij 269 enddo 270 ! write(std_out,*) itypat 271 ! do m=1,i 272 ! write(std_out,*) i,m 273 ! write(std_out,*) xij(i,m) 274 ! write(std_out,*) ij1,ij2 275 ! enddo 276 if(i==1) write(message,'(i3,14f7.3)') & 277 & i,xtemp, (hu(itypat)%uqmc(m),m=ij1,ij2) 278 if(i/=2*ndim.and.i/=1) write(message,'(i3,14f7.3)') i, & 279 & (hu(itypat)%uqmc(xij(i,m)), m=1,i-1),xtemp, (hu(itypat)%uqmc(m),m=ij1,ij2) 280 if(i==2*ndim) write(message,'(i3,14f7.3)') i, & 281 & (hu(itypat)%uqmc(xij(i,m)), m=1,i-1),xtemp 282 call wrtout(std_out, message,'COLL') 283 enddo 284 write(message,'(5x,a)') "--------------------------------------------------------" 285 call wrtout(std_out, message,'COLL') 286 ABI_FREE(xij) 287 else 288 hu(itypat)%upawu=zero 289 hu(itypat)%jpawu=zero 290 ! allocate(hu(itypat)%vee(0,0,0,0)) 291 endif 292 enddo ! itypat 293 294 end subroutine init_hu
m_hu/print_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
print_hu
FUNCTION
print density density interaction (used for DFT+DMFT)
INPUTS
ntypat = number of species prtopt = option for printing hu <type(hu_type)> = data for the interaction in DMFT.
OUTPUT
SOURCE
415 subroutine print_hu(hu,ntypat,prtopt) 416 417 use defs_basis 418 use m_crystal, only : crystal_t 419 implicit none 420 421 !Arguments ------------------------------------ 422 !type 423 integer, intent(in):: ntypat 424 type(hu_type),intent(in) :: hu(ntypat) 425 integer :: prtopt 426 427 !Local variables------------------------------- 428 integer :: itypat 429 integer :: lpawu,ms,ms1,m,ndim 430 character(len=500) :: message 431 ! ********************************************************************* 432 433 do itypat = 1 , ntypat 434 lpawu=hu(itypat)%lpawu 435 if(lpawu/=-1) then 436 ndim=2*lpawu+1 437 write(message,'(2a,i4)') ch10,' -------> For Correlated species' 438 call wrtout(std_out, message,'COLL') 439 if(prtopt==0) then 440 write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in Slm basis " 441 else if(prtopt==1) then 442 write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in diagonal basis" 443 else if(prtopt==2) then 444 write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in Ylm basis" 445 else if(prtopt==3) then 446 write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in JMJ basis" 447 endif 448 call wrtout(std_out, message,'COLL') 449 write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim) 450 call wrtout(std_out, message,'COLL') 451 do ms=1,2*ndim 452 write(message,'(i3,14f7.3)') & 453 & ms, (hu(itypat)%udens(ms,ms1),ms1=1,2*ndim) 454 call wrtout(std_out, message,'COLL') 455 enddo 456 write(message,'(5x,a)') "--------------------------------------------------------" 457 call wrtout(std_out, message,'COLL') 458 endif ! lpawu/=1 459 enddo ! ntypat 460 461 462 end subroutine print_hu
m_hu/printvee_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
printvee_hu
FUNCTION
print vee
INPUTS
vee = number of species lpawu = value of l
OUTPUT
SOURCE
1046 subroutine printvee_hu(ndim,vee,prtopt,basis,upawu,f2) 1047 1048 use defs_basis 1049 use m_crystal, only : crystal_t 1050 implicit none 1051 1052 !Arguments ------------------------------------ 1053 !type 1054 integer,intent(in) :: ndim,prtopt 1055 real(dp), intent(in) :: vee(ndim,ndim,ndim,ndim) 1056 real(dp), optional, intent(in) :: upawu,f2 1057 character(len=*), intent(in) :: basis 1058 1059 !Local variables------------------------------- 1060 integer :: abcomp,m1,m2,mi,mj,mk,ml 1061 real(dp),allocatable :: b0(:,:) 1062 real(dp),allocatable :: a2pp(:,:) 1063 real(dp),allocatable :: b2pp(:,:) 1064 real(dp):: aver52,aver72,avernondiag,averall,averu,averj,averurestricted 1065 character(len=2000) :: message 1066 ! ********************************************************************* 1067 write(message,'(5a)') ch10,& 1068 & ' Coulomb interaction in the ', trim(basis),' basis' 1069 call wrtout(std_out,message,'COLL') 1070 1071 if(prtopt>=2) then 1072 1073 write(message,'(2a)') ch10," <mi,mi|vee|mi mi> : U1" 1074 call wrtout(std_out, message,'COLL') 1075 do mi=1,ndim 1076 write(message,'(4i4,3x,e10.3)') mi,mi,mi,mi,vee(mi,mi,mi,mi) 1077 call wrtout(std_out, message,'COLL') 1078 enddo 1079 1080 write(message,'(2a)') ch10," <mi,mj|vee|mi mj> : U2" 1081 call wrtout(std_out, message,'COLL') 1082 do mi=1,ndim 1083 do mj=mi+1,ndim 1084 write(message,'(4i4,3x,e10.3)') mi,mj,mi,mj,vee(mi,mj,mi,mj) 1085 call wrtout(std_out, message,'COLL') 1086 enddo 1087 enddo 1088 1089 write(message,'(2a)') ch10," <mi,mj|vee|mj mi> : J" 1090 call wrtout(std_out, message,'COLL') 1091 do mi=1,ndim 1092 do mj=mi+1,ndim 1093 write(message,'(4i4,3x,e10.3)') mi,mj,mj,mi,vee(mi,mj,mj,mi) 1094 call wrtout(std_out, message,'COLL') 1095 enddo 1096 enddo 1097 1098 write(message,'(2a)') ch10," <mi,mi|vee|mj mj> : J" 1099 call wrtout(std_out, message,'COLL') 1100 do mi=1,ndim 1101 do mj=mi+1,ndim 1102 write(message,'(4i4,3x,e10.3)') mi,mi,mj,mj,vee(mi,mi,mj,mj) 1103 call wrtout(std_out, message,'COLL') 1104 enddo 1105 enddo 1106 1107 write(message,'(2a)') ch10," vee is non zero also for" 1108 call wrtout(std_out, message,'COLL') 1109 do mi=1,ndim 1110 do mj=1,ndim 1111 do mk=1,ndim 1112 do ml=1,ndim 1113 if((.not.(mi==mk.and.mj==ml)).and.(.not.(mi==ml.and.mj==mk)).and.& 1114 & .not.(mi==mj.and.mk==ml)) then 1115 if(vee(mi,mj,mk,ml)>tol8) then 1116 write(message,'(4i4,3x,e10.3)') mi,mj,mk,ml,vee(mi,mj,mk,ml) 1117 call wrtout(std_out, message,'COLL') 1118 endif 1119 endif 1120 enddo 1121 enddo 1122 enddo 1123 enddo 1124 write(message,'(a)') ch10 1125 call wrtout(std_out, message,'COLL') 1126 1127 endif 1128 1129 if(prtopt>=1) then 1130 1131 write(message,'(2x,a,3x,14f10.4)') "Um1m2=Vee(m1,m2,m1,m2)" 1132 call wrtout(std_out, message,'COLL') 1133 write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim) 1134 call wrtout(std_out, message,'COLL') 1135 do m1=1,ndim 1136 write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m1,m2),m2=1,ndim) 1137 call wrtout(std_out, message,'COLL') 1138 enddo 1139 write(message,'(a)') ch10; call wrtout(std_out, message,'COLL') 1140 1141 write(message,'(2x,a,3x,14f10.4)') "Jm1m2=Vee(m1,m2,m2,m1)" 1142 call wrtout(std_out, message,'COLL') 1143 write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim) 1144 call wrtout(std_out, message,'COLL') 1145 do m1=1,ndim 1146 write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m2,m1),m2=1,ndim) 1147 call wrtout(std_out, message,'COLL') 1148 enddo 1149 write(message,'(a)') ch10; call wrtout(std_out, message,'COLL') 1150 1151 write(message,'(2x,a,3x,14f10.4)') "Udens(m1,m2)" 1152 call wrtout(std_out, message,'COLL') 1153 write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim) 1154 call wrtout(std_out, message,'COLL') 1155 do m1=1,ndim 1156 write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1),m2=1,ndim) 1157 call wrtout(std_out, message,'COLL') 1158 enddo 1159 write(message,'(a)') ch10; call wrtout(std_out, message,'COLL') 1160 1161 if(prtopt>=3) then 1162 if(ndim==7) then 1163 averu=zero 1164 do m1=1,7 1165 do m2=1,7 1166 averu=vee(m1,m2,m1,m2)+averu 1167 enddo 1168 enddo 1169 averu=averu/49.d0 1170 averurestricted=zero 1171 do m1=1,7 1172 do m2=1,7 1173 if(m1/=m2) averurestricted=vee(m1,m2,m1,m2)+averurestricted 1174 enddo 1175 enddo 1176 averurestricted=averurestricted/42.d0 1177 averj=zero 1178 do m1=1,7 1179 do m2=1,7 1180 if(m1/=m2) averj=vee(m1,m2,m2,m1)+averj 1181 enddo 1182 enddo 1183 averj=averj/42 1184 averall=zero 1185 do m1=1,7 1186 do m2=1,7 1187 if(m1/=m2) averall=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+averall 1188 enddo 1189 enddo 1190 averall=averall/42 1191 !write(6,*) "averages Ylm U, U restricted,J, U-J",averu,averurestricted,averj,averall 1192 1193 endif 1194 1195 if(ndim==14.and.trim(basis)=='jmj') then 1196 aver52=zero 1197 do m1=1,6 1198 do m2=1,6 1199 aver52=vee(m1,m2,m1,m2)+aver52 1200 enddo 1201 enddo 1202 aver52=aver52/36.d0 1203 aver72=zero 1204 do m1=7,14 1205 do m2=7,14 1206 aver72=vee(m1,m2,m1,m2)+aver72 1207 enddo 1208 enddo 1209 aver72=aver72/64.d0 1210 avernondiag=zero 1211 do m1=1,6 1212 do m2=7,14 1213 avernondiag=vee(m1,m2,m1,m2)+avernondiag 1214 enddo 1215 enddo 1216 avernondiag=avernondiag/48.d0 1217 averall=zero 1218 do m1=1,14 1219 do m2=1,14 1220 averall=vee(m1,m2,m1,m2)+averall 1221 enddo 1222 enddo 1223 averall=averall/196.d0 1224 !write(6,*) "U averages",aver52,aver72,avernondiag,averall 1225 1226 1227 1228 aver52=zero 1229 do m1=1,6 1230 do m2=1,6 1231 if(m1/=m2) aver52=vee(m1,m2,m2,m1)+aver52 1232 enddo 1233 enddo 1234 aver52=aver52/30.d0 1235 aver72=zero 1236 do m1=7,14 1237 do m2=7,14 1238 if(m1/=m2) aver72=vee(m1,m2,m2,m1)+aver72 1239 enddo 1240 enddo 1241 aver72=aver72/56.d0 1242 avernondiag=zero 1243 do m1=1,6 1244 do m2=7,14 1245 avernondiag=vee(m1,m2,m2,m1)+avernondiag 1246 enddo 1247 enddo 1248 avernondiag=avernondiag/48.d0 1249 averall=zero 1250 do m1=1,14 1251 do m2=1,14 1252 if(m1/=m2) averall=vee(m1,m2,m2,m1)+averall 1253 enddo 1254 enddo 1255 averall=averall/182.d0 1256 !write(6,*) "J averages",aver52,aver72,avernondiag,averall 1257 1258 1259 1260 1261 aver52=zero 1262 do m1=1,6 1263 do m2=1,6 1264 if(m1/=m2) aver52=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+aver52 1265 enddo 1266 enddo 1267 aver52=aver52/30.d0 1268 aver72=zero 1269 do m1=7,14 1270 do m2=7,14 1271 if(m1/=m2) aver72=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+aver72 1272 enddo 1273 enddo 1274 aver72=aver72/56.d0 1275 avernondiag=zero 1276 do m1=1,6 1277 do m2=7,14 1278 avernondiag=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+avernondiag 1279 enddo 1280 enddo 1281 avernondiag=avernondiag/48.d0 1282 averall=zero 1283 do m1=1,14 1284 do m2=1,14 1285 if(m1/=m2) averall=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+averall 1286 enddo 1287 enddo 1288 averall=averall/182.d0 1289 !write(6,*) "U-J averages",aver52,aver72,avernondiag,averall 1290 1291 endif 1292 endif 1293 1294 if (present(upawu)) then 1295 ABI_MALLOC(a2pp,(ndim,ndim)) 1296 ABI_MALLOC(b2pp,(ndim,ndim)) 1297 ABI_MALLOC(b0,(ndim,ndim)) 1298 1299 ! write(message,'(2x,a,3x,14f10.4)') "For check with respect to Slater's paper" 1300 ! call wrtout(std_out, message,'COLL') 1301 ! ABI_MALLOC(f0,(ndim,ndim,ndim,ndim)) 1302 ! write(message,'(2x,a,3x,14f10.4)') "Vee(m1,m2,m1,m2)-F0*ao(m1,m2)" 1303 ! call wrtout(std_out, message,'COLL') 1304 ! write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim) 1305 ! call wrtout(std_out, message,'COLL') 1306 ! 1307 ! do m1=1,ndim 1308 ! write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m1,m2)-upawu,m2=1,ndim) 1309 ! call wrtout(std_out, message,'COLL') 1310 ! enddo 1311 ! 1312 ! f0=zero 1313 ! do m1=1,ndim 1314 ! f0(m1,m1,m1,m1)=upawu 1315 ! enddo 1316 ! write(message,'(a)') ch10 1317 ! call wrtout(std_out, message,'COLL') 1318 ! write(message,'(2x,a,3x,14f10.4)') "Vee(m1,m2,m2,m1)-F0*b0(m1,m2)" 1319 ! call wrtout(std_out, message,'COLL') 1320 ! write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim) 1321 ! call wrtout(std_out, message,'COLL') 1322 ! do m1=1,ndim 1323 ! write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m2,m1)-f0(m1,m2,m2,m1),m2=1,ndim) 1324 ! call wrtout(std_out, message,'COLL') 1325 ! enddo 1326 ! ABI_FREE(f0) 1327 1328 1329 b0=zero 1330 do m1=1,ndim 1331 b0(m1,m1)=upawu 1332 enddo 1333 abcomp=0 1334 if(ndim==3.and.present(f2).and.(trim(basis)=='slm')) then 1335 a2pp(:,:) = RESHAPE((/1,-2,1,-2,4,-2,-1,-2,-1/),(/3,3/)) 1336 a2pp=a2pp/25*f2+upawu 1337 b2pp(:,:) = RESHAPE((/1,3,6,3,4,3,6,3,1/),(/3,3/)) 1338 b2pp=b2pp/25*f2+b0 1339 abcomp=1 1340 else if(ndim==6.and.present(f2).and.(trim(basis)=='jmj')) then 1341 ABI_MALLOC(a2pp,(6,6)) 1342 ABI_MALLOC(b2pp,(6,6)) 1343 a2pp(:,:)=RESHAPE((/0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,-1,-1,1,0,& 1344 & 0,-1,1,1,-1,0,0,-1,1,1,-1,0,0,1,-1,-1,1/),(/6,6/)) 1345 a2pp=a2pp/25*f2+upawu 1346 b2pp(:,:)=RESHAPE((/0,0,1,2,3,4,0,0,4,3,2,1,1,4,1,2,2,0,2,3,& 1347 & 2,1,0,2,3,2,2,0,1,2,4,1,0,2,2,1/),(/6,6/)) 1348 b2pp=b2pp/25*f2+b0 1349 abcomp=1 1350 endif 1351 if(mod(ndim,3)==0.and.present(f2).and.abcomp==1) then 1352 write(message,'(2x,a)') "Exact result for Umm is" 1353 call wrtout(std_out,message,'COLL') 1354 do m1=1,ndim 1355 write(message,'(2x,i4,3x,14f10.4)') m1,(a2pp(m1,m2),m2=1,ndim) 1356 call wrtout(std_out, message,'COLL') 1357 enddo 1358 write(message,'(a)') ch10; call wrtout(std_out, message,'COLL') 1359 write(message,'(2x,a,3x,14f10.4)') "Exact result for Jmm is" 1360 call wrtout(std_out,message,'COLL') 1361 do m1=1,ndim 1362 write(message,'(2x,i4,3x,14f10.4)') m1,(b2pp(m1,m2),m2=1,ndim) 1363 call wrtout(std_out,message,'COLL') 1364 enddo 1365 write(message,'(a)') ch10; call wrtout(std_out, message,'COLL') 1366 endif 1367 ABI_FREE(a2pp) 1368 ABI_FREE(b2pp) 1369 ABI_FREE(b0) 1370 endif 1371 1372 endif 1373 1374 end subroutine printvee_hu
m_hu/reddd [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
reddd
FUNCTION
INPUTS
OUTPUT
SOURCE
1482 function reddd(mi,ndim) 1483 1484 use defs_basis 1485 implicit none 1486 1487 !Arguments ------------------------------------ 1488 !scalars 1489 integer,intent(in) :: mi,ndim 1490 integer :: reddd 1491 ! ************************************************************************* 1492 1493 if(mi<ndim+1) reddd=mi 1494 if(mi>=ndim+1) reddd=mi-ndim 1495 1496 end function reddd
m_hu/rotatevee_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
rotatevee_hu
FUNCTION
interaction udens in recomputed from new vee.
INPUTS
ntypat = number of species cryst_struc <type(crystal_t)>=crystal structure data
OUTPUT
SIDE EFFECT hu <type(hu_type)> = data for the interaction in DMFT.
SOURCE
567 subroutine rotatevee_hu(cryst_struc,hu,nspinor,nsppol,pawprtvol,rot_mat,udens_atoms,rot_type) 568 569 use defs_basis 570 use m_errors 571 572 use m_crystal, only : crystal_t 573 implicit none 574 575 !Arguments ------------------------------------ 576 !type 577 type(crystal_t),intent(in) :: cryst_struc 578 integer, intent(in):: nsppol,nspinor,pawprtvol,rot_type 579 type(hu_type),intent(in) :: hu(cryst_struc%ntypat) 580 type(coeff2c_type),optional,intent(in) :: rot_mat(cryst_struc%natom,nsppol) 581 type(coeff2_type),intent(inout) :: udens_atoms(cryst_struc%natom) 582 583 !Local variables------------------------------- 584 integer :: dim_vee,iatom,itypat 585 integer :: lpawu,m1,m2,m3,m4,mi,mj,mk,ml,natom,ndim,tndim 586 integer :: ms,ms1,nat_correl 587 character(len=500) :: message 588 character(len=30) :: basis_vee 589 real(dp) :: xsum,xsum2,xsumnew,xsum2new,uaver 590 complex(dpc),allocatable :: temp_mat(:,:) 591 complex(dpc),allocatable :: temp_mat2(:,:) 592 real(dp),allocatable :: veetemp(:,:,:,:),fk(:) 593 complex(dpc),allocatable :: veeylm(:,:,:,:),veeslm(:,:,:,:) 594 complex(dpc),allocatable :: veetmp(:,:,:),veetmp_slm(:,:,:),veetmp_ylm(:,:,:) 595 complex(dpc),allocatable :: veejmj(:,:,:,:),veerotated(:,:,:,:),veenew(:,:,:,:) 596 complex(dpc),allocatable :: veeylm2(:,:,:,:),veeslm2(:,:,:,:) 597 ! ********************************************************************* 598 natom=cryst_struc%natom 599 600 !================================================ 601 ! NSPINOR = 2 and J=0 602 !================================================ 603 if(hu(1)%jpawu_zero.and.nspinor==2) then 604 ! if(3==4) then 605 606 ! call vee2udens_hu(hu,cryst_struc%ntypat,2) 607 do iatom=1,natom 608 itypat=cryst_struc%typat(iatom) 609 lpawu=hu(itypat)%lpawu 610 if(lpawu.ne.-1) then 611 ndim=2*lpawu+1 612 write(message,'(2a,i4)') ch10,' -------> For Correlated Species', itypat 613 call wrtout(std_out, message,'COLL') 614 ! write(std_out,*)"ndim",ndim 615 ! write(std_out,'(a,20e10.4)')"udens before vee2udensaomt", udens_atoms(1)%value 616 ! write(std_out,*)"vee",hu(cryst_struc%typat(iatom))%vee 617 call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,hu(cryst_struc%typat(iatom))%vee,"Slm") 618 ! write(std_out,*)"udens after vee2udensatom", udens_atoms(1)%value 619 endif 620 enddo 621 ! write(std_out,*)"udensafter after", udens_atoms(1)%value 622 623 624 !================================================ 625 ! NSPINOR = 2 and J/=0 626 !================================================ 627 else if (nspinor==2) then 628 write(std_out,*) 629 write(std_out,*)" == Rotate interaction in the level basis " 630 631 632 !rot_type=0 ! keep original Slm basis 633 !rot_type=2 ! rotation to the Ylm basis ! for tests 634 !rot_type=4 ! use a rotation from Ylm basis to a new basis 635 !rot_type=3 ! rotation to the JmJ Basis ! for tests 636 !rot_type=1 ! use a rotation for diago of dmat,green, levels.. 637 638 do iatom=1,natom 639 itypat=cryst_struc%typat(iatom) 640 lpawu=hu(itypat)%lpawu 641 if(lpawu.ne.-1) then 642 ndim=2*lpawu+1 643 if(pawprtvol>=3) then 644 ! write(message,'(2a)') ch10," VEE INPUT AVANT TRANSFORMATION" 645 ! call wrtout(std_out, message,'COLL') 646 ! call printvee_hu(ndim,hu(itypat)%vee,1,'Slm') 647 endif 648 write(message,'(2a,i4)') ch10,' -------> For Correlated atom', iatom 649 call wrtout(std_out, message,'COLL') 650 tndim=nspinor*ndim 651 ABI_MALLOC(veejmj,(tndim,tndim,tndim,tndim)) 652 ABI_MALLOC(veeslm,(ndim,ndim,ndim,ndim)) 653 ABI_MALLOC(veetmp_slm,(ndim,ndim,4)) 654 ABI_MALLOC(veetmp_ylm,(ndim,ndim,4)) 655 ABI_MALLOC(veetmp,(ndim,ndim,4)) 656 ABI_MALLOC(veeslm2,(tndim,tndim,tndim,tndim)) 657 ABI_MALLOC(veeylm,(ndim,ndim,ndim,ndim)) 658 ABI_MALLOC(veeylm2,(tndim,tndim,tndim,tndim)) 659 ABI_MALLOC(veerotated,(tndim,tndim,tndim,tndim)) 660 ABI_MALLOC(fk,(0:lpawu)) 661 fk(0)=hu(itypat)%upawu 662 if (lpawu==1) then 663 fk(1)=hu(itypat)%jpawu*5 664 else if (lpawu==2) then 665 fk(1)=hu(itypat)%jpawu*14._dp/(One+hu(itypat)%f4of2_sla) 666 fk(2)=fk(1)*hu(itypat)%f4of2_sla 667 else if (lpawu==3) then 668 fk(1)=hu(itypat)%jpawu*6435._dp/(286._dp+195._dp*hu(itypat)%f4of2_sla+250._dp*hu(itypat)%f6of2_sla) 669 fk(2)=fk(1)*hu(itypat)%f4of2_sla 670 fk(3)=fk(1)*hu(itypat)%f6of2_sla 671 endif 672 673 674 675 676 ! ================================== 677 ! First define veeslm 678 ! ================================== 679 680 681 !! veeslm2(s1m1,s2m2,s3m3,s4m4)=vee(m1,m2,m3,m4)*delta_s1s3*delta_s2s4 682 call vee_ndim2tndim_hu(lpawu,hu(itypat)%vee,veeslm2,1) 683 684 !! veeslm(m1,m2,m3,m4)=cmplx(vee(m1,m2,m3,m4),zero) 685 ! ABI_FREE(veeslm) 686 veeslm(:,:,:,:)=cmplx(hu(itypat)%vee(:,:,:,:),zero) 687 !ABI_FREE(veeslm)! ici cela plante 688 689 !! build udens in the Slm basis and print it 690 call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,real(veeslm),"slm") 691 !ABI_FREE(veeslm) ! ici cela plante 692 693 dim_vee=ndim 694 basis_vee='Slm' 695 ! first print veeslm 696 !call printvee_hu(ndim,real(veeslm),1,basis_vee) 697 call printvee_hu(tndim,real(veeslm2),1,basis_vee) 698 699 ! ABI_FREE(veeslm) 700 701 ! ================================== 702 ! Then compute veerotated 703 ! ================================== 704 705 ! In the basis where levels/densitymatrix/green function is diagonalized 706 ! ================================================================================ 707 if (rot_type==1) then 708 ! --------------------- 709 710 veerotated=czero 711 do m1=1,tndim 712 do m2=1,tndim 713 do m3=1,tndim 714 do m4=1,tndim 715 do mi=1,tndim 716 do mj=1,tndim 717 do mk=1,tndim 718 do ml=1,tndim 719 veerotated(m1,m2,m3,m4)= veerotated(m1,m2,m3,m4) + & 720 & conjg(rot_mat(iatom,1)%value(mi,m1))* & 721 & conjg(rot_mat(iatom,1)%value(mj,m2))* & 722 & rot_mat(iatom,1)%value(mk,m3)* & 723 & rot_mat(iatom,1)%value(ml,m4)* & 724 & veeslm2(mi,mj,mk,ml) 725 enddo 726 enddo 727 enddo 728 enddo 729 enddo 730 enddo 731 enddo 732 enddo 733 734 dim_vee=tndim 735 basis_vee='Diagonal basis from Slm basis' 736 ! In the Ylm basis 737 ! ================================================================================ 738 else if (rot_type==2.or.rot_type==3.or.rot_type==4) then 739 ! --------------------------- 740 741 742 ! Change basis from slm to ylm basis 743 call vee_slm2ylm_hu(lpawu,veeslm,veeylm,1,2) 744 745 basis_vee='Ylm' 746 dim_vee=ndim 747 748 ! print interaction matrix in the ylm basis 749 call printvee_hu(ndim,real(veeylm),1,basis_vee,hu(itypat)%upawu) 750 751 ! print interaction matrix in the ylm basis from Slater tables 752 if(pawprtvol>=3) call udens_slatercondon_hu(fk,lpawu) 753 754 if (rot_type==3.or.rot_type==4) then 755 ! --------------------------- 756 ! 757 758 ! built large matrix 759 call vee_ndim2tndim_hu(lpawu,real(veeylm),veeylm2,1) 760 761 762 call printvee_hu(tndim,real(veeylm2),1,basis_vee) 763 764 endif 765 766 767 ! In the JmJ basis 768 ! ================================================================================ 769 if (rot_type==3) then 770 771 ! apply change of basis 772 call vee_ylm2jmj_hu(lpawu,veeylm2,veejmj,1) 773 774 ! print interaction matrix in the JMJ basis from Inglis and Julien tables 775 if(pawprtvol>=3) call udens_inglis_hu(fk,lpawu) 776 777 ! new dimension 778 779 780 781 dim_vee=tndim 782 basis_vee='JmJ' 783 784 else if (rot_type==4) then 785 786 call udens_inglis_hu(fk,lpawu) 787 veerotated=czero 788 789 do m1=1,tndim 790 do m2=1,tndim 791 do m3=1,tndim 792 do m4=1,tndim 793 do mi=1,tndim 794 do mj=1,tndim 795 do mk=1,tndim 796 do ml=1,tndim 797 veerotated(m1,m2,m3,m4)= veerotated(m1,m2,m3,m4) + & 798 & conjg(rot_mat(iatom,1)%value(mi,m1))* & 799 & conjg(rot_mat(iatom,1)%value(mj,m2))* & 800 & rot_mat(iatom,1)%value(mk,m3)* & 801 & rot_mat(iatom,1)%value(ml,m4)* & 802 & veeylm2(mi,mj,mk,ml) 803 enddo 804 enddo 805 enddo 806 enddo 807 enddo 808 enddo 809 enddo 810 enddo 811 812 dim_vee=tndim 813 basis_vee='Diagonal basis from Ylm' 814 815 endif 816 endif 817 ABI_MALLOC(veenew,(dim_vee,dim_vee,dim_vee,dim_vee)) 818 if(rot_type==0) then ; veenew=veeslm 819 else if(rot_type==1) then ; veenew=veerotated 820 else if(rot_type==2) then ; veenew=veeylm 821 else if(rot_type==3) then ; veenew=veejmj 822 else if(rot_type==4) then ; veenew=veerotated 823 endif 824 825 call printvee_hu(dim_vee,real(veenew),1,basis_vee,hu(itypat)%upawu,hu(itypat)%f2_sla) 826 ! call printvee_hu(dim_vee,real(veeylm),1,hu(itypat)%upawu) 827 828 uaver=zero 829 do ms=1,2*ndim 830 do ms1=1,2*ndim 831 udens_atoms(iatom)%value(ms,ms1)=veenew(ms,ms1,ms,ms1)-veenew(ms,ms1,ms1,ms) 832 uaver=uaver+udens_atoms(iatom)%value(ms,ms1) 833 enddo 834 enddo 835 836 call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,real(veenew),"basis_vee",prtonly=1) 837 838 839 ABI_FREE(veenew) 840 ABI_FREE(veeslm) 841 ABI_FREE(veeslm2) 842 ABI_FREE(veeylm) 843 ABI_FREE(veeylm2) 844 ABI_FREE(veejmj) 845 ABI_FREE(veerotated) 846 ABI_FREE(veetmp_slm) 847 ABI_FREE(veetmp) 848 ABI_FREE(veetmp_ylm) 849 ABI_FREE(fk) 850 endif 851 enddo 852 !ABI_ERROR("Aborting now!") 853 854 !================================================ 855 ! NSPINOR = 1 856 !================================================ 857 else if (nspinor==1) then 858 859 nat_correl=0 860 do iatom=1,natom 861 itypat=cryst_struc%typat(iatom) 862 lpawu=hu(itypat)%lpawu 863 if(lpawu.ne.-1) then 864 nat_correl=nat_correl+1 865 if(nat_correl>1.and.(hu(itypat)%jpawu>tol4)) then 866 write(message,'(3a)') ch10,' -------> Warning: several atoms: '& 867 & ,' not extensively tested ' 868 ABI_WARNING(message) 869 endif 870 871 ! ! ================================================================ 872 ! ! If rotation for spin 2 and rotation for spin 1 are not equal print 873 ! ! then print a warning 874 ! ! useful only for magnetic case 875 ! ! ================================================================ 876 ndim=2*lpawu+1 877 tndim=nspinor*ndim 878 do m1=1,tndim 879 do m2=1,tndim 880 if(nsppol==2) then 881 if(abs(rot_mat(iatom,1)%value(m1,m2)-rot_mat(iatom,2)%value(m1,m2))>tol4.and.pawprtvol>=3) then 882 write(message,'(2a,i4)') ch10,' rot_mat differs for value of isppol but value for isppol=2 not used' 883 call wrtout(std_out, message,'COLL') 884 write(message,'(a,4e16.8)') ch10,rot_mat(iatom,1)%value(m1,m2),rot_mat(iatom,2)%value(m1,m2) 885 call wrtout(std_out, message,'COLL', do_flush=.True.) 886 endif 887 endif 888 end do 889 end do 890 ABI_MALLOC(temp_mat,(ndim,ndim)) 891 ABI_MALLOC(temp_mat2,(ndim,ndim)) 892 ABI_MALLOC(veetemp,(ndim,ndim,ndim,ndim)) 893 temp_mat(:,:)=czero 894 temp_mat2(:,:)=czero 895 896 ! ! ================================================= 897 ! ! See if rotation is complex or real 898 ! ! ================================================= 899 do mi=1,ndim 900 do m1=1,ndim 901 if(abs(aimag(rot_mat(iatom,1)%value(mi,m1)))>tol8.and.pawprtvol>=3) then 902 write(message,'(2a,2i6,2e14.3)') ch10,"rot_mat is complex for", & 903 & mi,m1,rot_mat(iatom,1)%value(mi,m1) 904 call wrtout(std_out, message,'COLL') 905 endif 906 enddo 907 enddo 908 909 910 911 ! ! write vee for information with a classification. 912 if(pawprtvol>=3) then 913 call printvee_hu(ndim,hu(itypat)%vee,2,'original') 914 endif 915 916 ! ! Compute rotated vee. 917 veetemp=zero 918 do m1=1,ndim 919 do m2=1,ndim 920 do m3=1,ndim 921 do m4=1,ndim 922 do mi=1,ndim 923 do mj=1,ndim 924 do mk=1,ndim 925 do ml=1,ndim 926 ! if((mi==mk.and.mj==ml).or.(mi==ml.and.mj==mk)) then 927 veetemp(m1,m2,m3,m4)= veetemp(m1,m2,m3,m4) + & 928 & real( & 929 & conjg(rot_mat(iatom,1)%value(mi,m1))* & 930 & conjg(rot_mat(iatom,1)%value(mj,m2))* & 931 & rot_mat(iatom,1)%value(mk,m3)* & 932 & rot_mat(iatom,1)%value(ml,m4)* & 933 & hu(itypat)%vee(mi,mj,mk,ml)& 934 ) 935 ! endif 936 enddo 937 enddo 938 enddo 939 enddo 940 enddo 941 enddo 942 enddo 943 enddo 944 ABI_FREE(temp_mat) 945 ABI_FREE(temp_mat2) 946 xsum=zero 947 xsum2=zero 948 xsumnew=zero 949 xsum2new=zero 950 do m1=1,ndim 951 do m2=1,ndim 952 xsum=xsum+hu(itypat)%vee(m1,m2,m1,m2) 953 xsum2=xsum2+hu(itypat)%vee(m1,m2,m2,m1) 954 xsumnew=xsumnew+veetemp(m1,m2,m1,m2) 955 xsum2new=xsum2new+veetemp(m1,m2,m2,m1) 956 enddo 957 enddo 958 if(abs(xsum-xsumnew)>tol5.or.abs(xsum2-xsum2new)>tol5) then 959 write(message,'(2a)') ch10," BUG: New interaction after rotation do not respect sum rules" 960 call wrtout(std_out, message,'COLL') 961 write(message,'(2a,2f14.3)') ch10,' Comparison of \sum_{m1,m3} vee(m1,m3,m1,m3) before and after rotation is',& 962 & xsum,xsumnew 963 call wrtout(std_out, message,'COLL') 964 write(message,'(2a,2f14.3)') ch10,' Comparison of \sum_{m1,m3} vee(m1,m3,m3,m1) before and after rotation is',& 965 & xsum2,xsum2new 966 call wrtout(std_out, message,'COLL') 967 endif 968 if(pawprtvol>=3) then 969 write(message,'(2a)') ch10," VEE ROTATED" 970 call wrtout(std_out, message,'COLL') 971 call printvee_hu(tndim,veetemp,2,'diag1') 972 write(message,'(a)') ch10 973 call wrtout(std_out, message,'COLL') 974 endif 975 976 write(message,'(2a,i4)') ch10,' -------> For Correlated Species', itypat 977 call wrtout(std_out, message,'COLL') 978 979 call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,veetemp,"diag") 980 981 ! udens_atoms(iatom)%value=zero 982 ! ij=0 983 ! do ms=1,2*ndim-1 984 ! do ms1=ms+1,2*ndim 985 ! ij=ij+1 986 ! if(ms<=ndim.and.ms1>ndim) then 987 ! m1 = ms1 - ndim 988 ! m = ms 989 ! hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1) 990 ! udens_atoms(iatom)%value(ms,ms1)= veetemp(m,m1,m,m1) 991 ! udens_atoms(iatom)%value(ms1,ms)= udens_atoms(iatom)%value(ms,ms1) 992 ! else if(ms<=ndim.and.ms1<=ndim) then 993 ! m1 = ms1 994 ! m = ms 995 ! hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m) 996 ! udens_atoms(iatom)%value(ms,ms1)= hu(itypat)%uqmc(ij) 997 ! udens_atoms(iatom)%value(ms1,ms)= udens_atoms(iatom)%value(ms,ms1) 998 ! else 999 ! m1 = ms1 - ndim 1000 ! m = ms - ndim 1001 ! hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m) 1002 ! udens_atoms(iatom)%value(ms,ms1)= hu(itypat)%uqmc(ij) 1003 ! udens_atoms(iatom)%value(ms1,ms)= udens_atoms(iatom)%value(ms,ms1) 1004 ! endif 1005 ! enddo 1006 ! enddo 1007 ! write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation " 1008 ! call wrtout(std_out, message,'COLL') 1009 ! write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim) 1010 ! call wrtout(std_out, message,'COLL') 1011 ! do ms=1,2*ndim 1012 ! write(message,'(i3,14f7.3)') & 1013 ! & ms, (udens_atoms(iatom)%value(ms,ms1),ms1=1,2*ndim) 1014 ! call wrtout(std_out, message,'COLL') 1015 ! enddo 1016 ! write(message,'(5x,a)') "--------------------------------------------------------" 1017 ! call wrtout(std_out, message,'COLL') 1018 ABI_FREE(veetemp) 1019 endif ! lpawu/=1 1020 ! call print_hu(hu,cryst_struc%ntypat,1) 1021 1022 enddo ! iatom 1023 ! call print_hu(hu,cryst_struc%ntypat,1) 1024 ! call vee2udens_hu(hu,cryst_struc%ntypat,2) 1025 endif ! nspinor==1 1026 1027 1028 end subroutine rotatevee_hu
m_hu/vee2udens_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
print_hu
FUNCTION
interaction udens in recomputed from new vee.
INPUTS
ntypat = number of species prtopt = option for printing
OUTPUT
SIDE EFFECT hu <type(hu_type)> = data for the interaction in DMFT.
SOURCE
483 subroutine vee2udens_hu(hu,ntypat,prtopt) 484 485 use defs_basis 486 use m_crystal, only : crystal_t 487 implicit none 488 489 !Arguments ------------------------------------ 490 !type 491 integer, intent(in):: ntypat 492 type(hu_type),intent(inout) :: hu(ntypat) 493 integer :: prtopt 494 495 !Local variables------------------------------- 496 integer :: ij,itypat 497 integer :: lpawu,m1,ms,ms1,m,ndim 498 character(len=500) :: message 499 ! ********************************************************************* 500 do itypat=1,ntypat 501 lpawu=hu(itypat)%lpawu 502 if(lpawu.ne.-1) then 503 ndim=2*lpawu+1 504 write(message,'(2a,i4)') ch10,' -------> For Correlated Species', itypat 505 call wrtout(std_out, message,'COLL') 506 507 hu(itypat)%udens=zero 508 ij=0 509 do ms=1,2*ndim-1 510 ! xij(ms,ms)=0 511 do ms1=ms+1,2*ndim 512 ij=ij+1 513 ! xij(ms,ms1)=ij 514 ! xij(ms1,ms)=ij 515 if(ms<=ndim.and.ms1>ndim) then 516 m1 = ms1 - ndim 517 m = ms 518 hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1) 519 hu(itypat)%udens(ms,ms1)= hu(itypat)%vee(m,m1,m,m1) 520 hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1) 521 else if(ms<=ndim.and.ms1<=ndim) then 522 m1 = ms1 523 m = ms 524 hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m) 525 hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij) 526 hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1) 527 else 528 m1 = ms1 - ndim 529 m = ms - ndim 530 hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m) 531 hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij) 532 hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1) 533 endif 534 enddo 535 enddo 536 ! xij(2*ndim,2*ndim)=0 537 ! write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation " 538 ! call wrtout(std_out, message,'COLL') 539 ! write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim) 540 ! call wrtout(std_out, message,'COLL') 541 endif 542 enddo ! itypat 543 call print_hu(hu,ntypat,prtopt) 544 545 546 end subroutine vee2udens_hu
m_hu/vee2udensatom_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
vee2udensatom_hu
FUNCTION
compute density density interaction (used for DFT+DMFT)
INPUTS
ntypat = number of species hu <type(hu_type)> = data for the interaction in DMFT.
OUTPUT
SOURCE
1392 subroutine vee2udensatom_hu(ndim,nspinor,udens_atoms,veetemp,basis,prtonly) 1393 1394 use defs_basis 1395 use m_crystal, only : crystal_t 1396 implicit none 1397 1398 !Arguments ------------------------------------ 1399 !type 1400 integer, intent(in):: ndim,nspinor 1401 real(dp),intent(out) :: udens_atoms(2*ndim,2*ndim) 1402 !real(dp), intent(in) :: veetemp(nspinor*ndim,nspinor*ndim,nspinor*ndim,nspinor*ndim) 1403 real(dp), intent(in) :: veetemp(ndim,ndim,ndim,ndim) 1404 character(len=*), intent(in) :: basis 1405 integer, intent(in), optional:: prtonly 1406 1407 !Local variables------------------------------- 1408 integer :: ij,ms,ms1,m,m1,prt 1409 character(len=1000) :: message 1410 ! ********************************************************************* 1411 prt=1 1412 if(present(prtonly)) prt=2 1413 if (nspinor==1.or.nspinor==2.and.prt<=1) then 1414 udens_atoms=zero 1415 ij=0 1416 do ms=1,2*ndim-1 1417 do ms1=ms+1,2*ndim 1418 ij=ij+1 1419 if(ms<=ndim.and.ms1>ndim) then 1420 m1 = ms1 - ndim 1421 m = ms 1422 ! hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1) 1423 udens_atoms(ms,ms1)= veetemp(m,m1,m,m1) 1424 udens_atoms(ms1,ms)= udens_atoms(ms,ms1) 1425 ! write(6,*)"A", ms,ms1,udens_atoms(ms,ms1) 1426 else if(ms<=ndim.and.ms1<=ndim) then 1427 m1 = ms1 1428 m = ms 1429 ! hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m) 1430 udens_atoms(ms,ms1)= veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m) 1431 udens_atoms(ms1,ms)= udens_atoms(ms,ms1) 1432 ! write(6,*)"B", ms,ms1,udens_atoms(ms,ms1) 1433 else 1434 m1 = ms1 - ndim 1435 m = ms - ndim 1436 ! hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m) 1437 udens_atoms(ms,ms1)= veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m) 1438 udens_atoms(ms1,ms)= udens_atoms(ms,ms1) 1439 ! write(6,*)"C", ms,ms1,udens_atoms(ms,ms1) 1440 endif 1441 enddo 1442 enddo 1443 1444 ! else if(nspinor==2) then 1445 ! 1446 ! do ms=1,2*ndim 1447 ! do ms1=1,2*ndim 1448 ! udens_atoms(ms,ms1)=veetemp(ms,ms1,ms,ms1)-veetemp(ms,ms1,ms1,ms) 1449 ! enddo 1450 ! enddo 1451 1452 endif 1453 1454 1455 write(message,'(4a)') ch10,"-------- Interactions in the ",trim(basis)," basis " 1456 call wrtout(std_out, message,'COLL') 1457 write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim) 1458 call wrtout(std_out, message,'COLL') 1459 do ms=1,2*ndim 1460 write(message,'(i3,14f7.3)') & 1461 & ms, (udens_atoms(ms,ms1),ms1=1,2*ndim) 1462 call wrtout(std_out, message,'COLL') 1463 enddo 1464 write(message,'(5x,a)') "--------------------------------------------------------" 1465 call wrtout(std_out, message,'COLL') 1466 1467 end subroutine vee2udensatom_hu
m_hu/vee_ndim2tndim_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
vee_ndim2tndim_hu
FUNCTION
Change a matrix of interaction of dimension [(2*lcor+1)]**2 into a full spin and orbital interaction matrix of dimension [2*(2l+1)]**4
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (BA) 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
lcor= angular momentum, size of the matrix is 2(2*lcor+1) mat_inp_c= real input matrix prtvol=printing volume option= 1 : Vout_(s1m1,s2m2,s3m3,s4m4)=Vinp_(m1,m2,m3,m4)*delta_s1s3*delta_s2s4
OUTPUT
mat_out_c= Complex output matrix
NOTES
SOURCE
1752 subroutine vee_ndim2tndim_hu(lcor,mat_inp_c,mat_out_c,option) 1753 1754 use defs_basis 1755 use m_errors 1756 use m_abicore 1757 implicit none 1758 1759 !Arguments --------------------------------------------- 1760 !scalars 1761 integer,intent(in) :: lcor,option 1762 !arrays 1763 real(dp), intent(in) :: mat_inp_c(2*lcor+1,2*lcor+1,2*lcor+1,2*lcor+1) 1764 complex(dpc), intent(out) :: mat_out_c(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1)) 1765 1766 !Local variables --------------------------------------- 1767 !scalars 1768 integer :: m1,m2,m3,m4,is1,is2,is3,is4,ndim,s1,s2,s3,s4 1769 1770 ! ********************************************************************* 1771 ndim=2*lcor+1 1772 mat_out_c=czero 1773 if(option==1) then 1774 do m1=1,ndim 1775 do m2=1,ndim 1776 do m3=1,ndim 1777 do m4=1,ndim 1778 do is1=1,2 1779 do is2=1,2 1780 1781 is3=is1 ; is4=is2 1782 1783 s1=(is1-1)*ndim ; s2=(is2-1)*ndim ; s3=(is3-1)*ndim ; s4=(is4-1)*ndim 1784 1785 mat_out_c(m1+s1,m2+s2,m3+s3,m4+s4)= cmplx(mat_inp_c(m1,m2,m3,m4),zero) 1786 1787 enddo 1788 enddo 1789 enddo 1790 enddo 1791 enddo 1792 enddo 1793 endif 1794 1795 1796 end subroutine vee_ndim2tndim_hu
m_hu/vee_ndim2tndim_hu_r [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
vee_ndim2tndim_hu_r
FUNCTION
Change a matrix of interaction of dimension [(2*lcor+1)]**2 into a full spin and orbital interaction matrix of dimension [2*(2l+1)]**4
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (BA) 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
lcor= angular momentum, size of the matrix is 2(2*lcor+1) mat_inp_c= real input matrix prtvol=printing volume option= 1 : Vout_(s1m1,s2m2,s3m3,s4m4)=Vinp_(m1,m2,m3,m4)*delta_s1s3*delta_s2s4
OUTPUT
mat_out_c= real output matrix
NOTES
SOURCE
1675 subroutine vee_ndim2tndim_hu_r(lcor,mat_inp_c,mat_out_c,option) 1676 1677 use defs_basis 1678 use m_errors 1679 use m_abicore 1680 implicit none 1681 1682 !Arguments --------------------------------------------- 1683 !scalars 1684 integer,intent(in) :: lcor,option 1685 !arrays 1686 real(dp), intent(in) :: mat_inp_c(:,:,:,:) !real(dp), intent(inout) :: mat_inp_c(:,:,:,:) !(2*lcor+1,2*lcor+1,2*lcor+1,2*lcor+1) real(dp), intent(in), pointer :: mat_inp_c(:,:,:,:) 1687 real(dp), intent(out) :: mat_out_c(:,:,:,:) !(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1)) 1688 1689 !Local variables --------------------------------------- 1690 !scalars 1691 integer :: m1,m2,m3,m4,is1,is2,is3,is4,ndim,s1,s2,s3,s4 1692 1693 ! ********************************************************************* 1694 ndim=2*lcor+1 1695 mat_out_c=czero 1696 1697 if(option==1) then 1698 do m1=1,ndim 1699 do m2=1,ndim 1700 do m3=1,ndim 1701 do m4=1,ndim 1702 do is1=1,2 1703 do is2=1,2 1704 1705 is3=is1 ; is4=is2 1706 1707 s1=(is1-1)*ndim ; s2=(is2-1)*ndim ; s3=(is3-1)*ndim ; s4=(is4-1)*ndim 1708 1709 mat_out_c(m1+s1,m2+s2,m3+s3,m4+s4)= mat_inp_c(m1,m2,m3,m4) 1710 1711 enddo 1712 enddo 1713 enddo 1714 enddo 1715 enddo 1716 enddo 1717 endif 1718 1719 1720 end subroutine vee_ndim2tndim_hu_r
m_hu/vee_slm2ylm_hu [ Functions ]
[ Top ] [ m_hu ] [ Functions ]
NAME
vee_slm2ylm_hu
FUNCTION
For a given angular momentum lcor, change a matrix of interaction of dimension (2*lcor+1) from the Slm to the Ylm basis if option==1 or from Ylm to Slm if !option==2
COPYRIGHT
Copyright (C) 1998-2022 ABINIT group (BA) 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
lcor= angular momentum, size of the matrix is 2(2*lcor+1) mat_inp_c= Input matrix option= 1 Change matrix from Slm to Ylm basis 2 Change matrix from Ylm to Slm basis optspin= 1 Spin up are first 2 Spin dn are first prtvol=printing volume
OUTPUT
mat_inp_c= Output matrix in Ylm or Slm basis according to option
NOTES
SOURCE
1530 subroutine vee_slm2ylm_hu(lcor,mat_inp_c,mat_out_c,option,prtvol) 1531 1532 use defs_basis 1533 use m_errors 1534 use m_abicore 1535 implicit none 1536 1537 !Arguments --------------------------------------------- 1538 !scalars 1539 integer,intent(in) :: lcor,option,prtvol 1540 !arrays 1541 complex(dpc), intent(in) :: mat_inp_c(2*lcor+1,2*lcor+1,2*lcor+1,2*lcor+1) 1542 complex(dpc), intent(out) :: mat_out_c(2*lcor+1,2*lcor+1,2*lcor+1,2*lcor+1) 1543 1544 !Local variables --------------------------------------- 1545 !scalars 1546 integer :: gm,hm,jm,gg,hh,ii,jj,ll,mm,im 1547 real(dp),parameter :: invsqrt2=one/sqrt2 1548 real(dp) :: onem 1549 complex(dpc) :: tmp2 1550 character(len=500) :: message 1551 !arrays 1552 complex(dpc),allocatable :: slm2ylm(:,:) 1553 1554 ! ********************************************************************* 1555 1556 1557 if (option/=1.and.option/=2) then 1558 message=' option=/1 or 2 !' 1559 ABI_BUG(message) 1560 end if 1561 1562 if(abs(prtvol)>2) then 1563 write(message,'(3a)') ch10, " vee_slm2ylm_hu" 1564 call wrtout(std_out,message,'COLL') 1565 end if 1566 1567 if(abs(prtvol)>2) then 1568 if(option==1) then 1569 write(message,'(3a)') ch10,"matrix in Slm basis is changed into Ylm basis" 1570 call wrtout(std_out,message,'COLL') 1571 else if(option==2) then 1572 write(message,'(3a)') ch10,"matrix in Ylm basis is changed into Slm basis" 1573 call wrtout(std_out,message,'COLL') 1574 end if 1575 end if 1576 1577 ll=lcor 1578 ABI_MALLOC(slm2ylm,(2*ll+1,2*ll+1)) 1579 slm2ylm=czero 1580 mat_out_c=czero 1581 1582 ! ===== Definitions of slm2ylm 1583 do im=1,2*ll+1 1584 mm=im-ll-1;jm=-mm+ll+1 ! mmj=-mm 1585 !! im is in {1,....2*ll+1} 1586 !! mm is in {-ll,....+ll} 1587 !! jm is in {2*ll+1,....,1} 1588 onem=dble((-1)**mm) 1589 if (mm> 0) then ! im in {ll+1,2ll+1} and jm in {ll+1,1} 1590 slm2ylm(im,im)= cmplx(onem*invsqrt2,zero,kind=dp) 1591 slm2ylm(jm,im)= cmplx(invsqrt2, zero,kind=dp) 1592 end if 1593 if (mm==0) then 1594 slm2ylm(im,im)=cone 1595 end if 1596 if (mm< 0) then 1597 slm2ylm(im,im)= cmplx(zero, invsqrt2,kind=dp) 1598 slm2ylm(jm,im)=-cmplx(zero,onem*invsqrt2,kind=dp) 1599 end if 1600 end do 1601 ! do im=1,2*ll+1 1602 ! write(message,'(7(2f14.5))') (slm2ylm(im,jm),jm=1,2*ll+1) 1603 ! call wrtout(std_out,message,'COLL') 1604 ! end do 1605 1606 ! ===== Definitions of slm2ylm 1607 !!!! pawtab(itypat)%vee(m11,m31,m21,m41)= <m11 m31| vee| m21 m41 > 1608 !!!! pawtab(itypat)%vee(m11,m21,m31,m41)= <m11 m21| vee| m31 m41 > 1609 1610 do jm=1,2*ll+1 1611 do im=1,2*ll+1 1612 do hm=1,2*ll+1 1613 do gm=1,2*ll+1 1614 tmp2=czero 1615 do gg=1,2*ll+1 1616 do hh=1,2*ll+1 1617 do ii=1,2*ll+1 1618 do jj=1,2*ll+1 1619 if(option==1) then 1620 tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*(slm2ylm(im,ii))*CONJG(slm2ylm(jm,jj))& 1621 & *(slm2ylm(gm,gg))*CONJG(slm2ylm(hm,hh)) 1622 ! if(gm==1.and.hm==1.and.im==1.and.jm==1) then 1623 ! write(6,'(4i4,2f10.5,2f10.5)') gg,hh,ii,jj,tmp2,mat_inp_c(gg,hh,ii,jj) 1624 ! write(6,*) "i1" 1625 ! endif 1626 else if(option==2) then 1627 tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*CONJG(slm2ylm(ii,im))*(slm2ylm(jj,jm))& 1628 & *CONJG(slm2ylm(gg,gm))*(slm2ylm(hh,hm)) 1629 end if 1630 end do 1631 end do 1632 end do 1633 end do 1634 ! mat_out_c(gm,hm,im,jm)=tmp2 1635 mat_out_c(gm,im,hm,jm)=tmp2 1636 end do 1637 end do 1638 end do 1639 end do 1640 1641 ABI_FREE(slm2ylm) 1642 1643 end subroutine vee_slm2ylm_hu