TABLE OF CONTENTS


ABINIT/m_hu [ Modules ]

[ Top ] [ Modules ]

NAME

  m_hu

FUNCTION

COPYRIGHT

 Copyright (C) 2006-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 .

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

24 #if defined HAVE_CONFIG_H
25 #include "config.h"
26 #endif
27 
28 
29 #include "abi_common.h"
30 
31 MODULE m_hu
32 
33  use m_abicore
34 
35  use defs_basis
36  use m_pawtab, only : pawtab_type
37 
38  implicit none
39 
40  private
41 
42  public :: init_hu
43  public :: destroy_hu
44 ! public :: qmc_hu
45  public :: print_hu
46  public :: vee2udens_hu
47  public :: rotatevee_hu
48  public :: printvee_hu
49  public :: vee2udensatom_hu
50  public :: vee_slm2ylm_hu
51  public :: vee_ndim2tndim_hu
52  public :: vee_ndim2tndim_hu_r
53  public :: udens_slatercondon_hu
54  public :: udens_inglis_hu

ABINIT/udens_inglis_hu [ Functions ]

[ Top ] [ 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-2018 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

2363 subroutine udens_inglis_hu(fk,lcor)
2364 
2365  use defs_basis
2366  use m_errors
2367  use m_abicore
2368 
2369 !This section has been created automatically by the script Abilint (TD).
2370 !Do not modify the following lines by hand.
2371 #undef ABI_FUNC
2372 #define ABI_FUNC 'udens_inglis_hu'
2373 !End of the abilint section
2374 
2375  implicit none
2376 
2377 !Arguments ---------------------------------------------
2378 !scalars
2379  integer,intent(in) :: lcor
2380  real(dp), intent(in) :: fk(0:lcor)
2381 
2382 !Local variables ---------------------------------------
2383 !scalars
2384  character(len=500) :: message
2385  integer :: m1,m2,kk,tndim
2386 !arrays
2387  real(dp), allocatable :: udens(:,:)
2388  real(dp), allocatable :: jdens(:,:)
2389  real(dp),allocatable :: app(:,:,:)
2390  real(dp),allocatable :: bpp(:,:,:)
2391  real(dp),allocatable :: a2pp(:,:)
2392  real(dp),allocatable :: b2pp(:,:)
2393 
2394 !*********************************************************************
2395  tndim=2*(2*lcor+1)
2396  ABI_ALLOCATE(app,(0:lcor,tndim,tndim))
2397  ABI_ALLOCATE(bpp,(0:lcor,tndim,tndim))
2398  ABI_ALLOCATE(a2pp,(tndim,tndim))
2399  ABI_ALLOCATE(b2pp,(tndim,tndim))
2400 
2401  ABI_ALLOCATE(udens,(tndim,tndim))
2402  ABI_ALLOCATE(jdens,(tndim,tndim))
2403  udens=zero
2404  jdens=zero
2405  a2pp=zero
2406  b2pp=zero
2407  app=zero
2408  bpp=zero
2409  if(lcor==1) then
2410    app(0,:,:)=one
2411    a2pp(:,:)=RESHAPE((/0._dp,0._dp, 0._dp, 0._dp, 0._dp, 0._dp,&
2412 &                      0._dp,0._dp, 0._dp, 0._dp, 0._dp, 0._dp,&
2413 &                      0._dp,0._dp, 1._dp,-1._dp,-1._dp, 1._dp,&
2414 &                      0._dp,0._dp,-1._dp, 1._dp, 1._dp,-1._dp,&
2415 &                      0._dp,0._dp,-1._dp, 1._dp, 1._dp,-1._dp,&
2416 &                      0._dp,0._dp, 1._dp,-1._dp,-1._dp, 1._dp/),(/6,6/))
2417    app(1,:,:)=a2pp(:,:)
2418    app(1,:,:)=app(1,:,:)/25._dp
2419 
2420    b2pp(:,:)=RESHAPE((/1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,&
2421 &                      0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,&
2422 &                      0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,&
2423 &                      0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,&
2424 &                      0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,&
2425 &                      0._dp,0._dp,0._dp,0._dp,0._dp,1._dp /),(/6,6/))
2426    bpp(0,:,:)=b2pp(:,:)
2427    b2pp(:,:)=RESHAPE((/0._dp,0._dp,1._dp,2._dp,3._dp,4._dp,&
2428 &                      0._dp,0._dp,4._dp,3._dp,2._dp,1._dp,&
2429 &                      1._dp,4._dp,1._dp,2._dp,2._dp,0._dp,&
2430 &                      2._dp,3._dp,2._dp,1._dp,0._dp,2._dp,&
2431 &                      3._dp,2._dp,2._dp,0._dp,1._dp,2._dp,&
2432 &                      4._dp,1._dp,0._dp,2._dp,2._dp,1._dp /),(/6,6/))
2433    bpp(1,:,:)=b2pp(:,:)
2434    bpp(1,:,:)=bpp(1,:,:)/25._dp
2435  endif
2436  if(lcor==2) then
2437    app(0,:,:)=one
2438    a2pp(:,:)=RESHAPE((/ 49._dp,-49._dp,-49._dp, 49._dp, 70._dp,-14._dp,-56._dp,-56._dp,-14._dp, 70._dp,&
2439 &                      -49._dp, 49._dp, 49._dp,-49._dp,-70._dp, 14._dp, 56._dp, 56._dp, 14._dp,-70._dp,&
2440 &                      -49._dp, 49._dp, 49._dp,-49._dp,-70._dp, 14._dp, 56._dp, 56._dp, 14._dp,-70._dp,&
2441 &                       49._dp,-49._dp,-49._dp, 49._dp, 70._dp,-14._dp,-56._dp,-56._dp,-14._dp, 70._dp,&
2442 &                       70._dp,-70._dp,-70._dp, 70._dp,100._dp,-20._dp,-80._dp,-80._dp,-20._dp,100._dp,&
2443 &                      -14._dp, 14._dp, 14._dp,-14._dp,-20._dp,  4._dp, 16._dp, 16._dp,  4._dp,-20._dp,&
2444 &                      -56._dp, 56._dp, 56._dp,-56._dp,-80._dp, 16._dp, 64._dp, 64._dp, 16._dp,-80._dp,&
2445 &                      -56._dp, 56._dp, 56._dp,-56._dp,-80._dp, 16._dp, 64._dp, 64._dp, 16._dp,-80._dp,&
2446 &                      -14._dp, 14._dp, 14._dp,-14._dp,-20._dp,  4._dp, 16._dp, 16._dp,  4._dp,-20._dp,&
2447 &                       70._dp,-70._dp,-70._dp, 70._dp,100._dp,-20._dp,-80._dp,-80._dp,-20._dp,100._dp/),(/10,10/))
2448    app(1,:,:)=a2pp(:,:)/1225._dp
2449    app(2,:,:)=zero
2450 
2451    b2pp(:,:)=RESHAPE((/1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,&
2452 &                      0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,&
2453 &                      0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,&
2454 &                      0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,&
2455 &                      0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,0._dp,&
2456 &                      0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,0._dp,&
2457 &                      0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,0._dp,&
2458 &                      0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,0._dp,&
2459 &                      0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp,0._dp,&
2460 &                      0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,0._dp,1._dp/),(/10,10/))
2461    bpp(0,:,:)=b2pp(:,:)
2462    b2pp(:,:)=RESHAPE((/ 49._dp, 98._dp, 98._dp,  0._dp, 30._dp, 36._dp, 27._dp, 12._dp,  0._dp,  0._dp,&
2463 &                       98._dp, 49._dp,  0._dp, 98._dp, 40._dp,  2._dp,  6._dp, 25._dp, 32._dp,  0._dp,&
2464 &                       98._dp,  0._dp, 49._dp, 98._dp,  0._dp, 32._dp, 25._dp,  6._dp,  2._dp, 40._dp,&
2465 &                        0._dp, 98._dp, 98._dp, 49._dp,  0._dp,  0._dp, 12._dp, 27._dp, 36._dp, 30._dp,&
2466 &                       30._dp, 40._dp,  0._dp,  0._dp,100._dp,120._dp, 60._dp,  0._dp,  0._dp,  0._dp,&
2467 &                       36._dp,  2._dp, 32._dp,  0._dp,120._dp,  4._dp, 48._dp,108._dp,  0._dp,  0._dp,&
2468 &                       27._dp,  6._dp, 25._dp, 12._dp, 60._dp, 48._dp, 64._dp,  0._dp,108._dp,  0._dp,&
2469 &                       12._dp, 25._dp,  6._dp, 27._dp,  0._dp,108._dp,  0._dp, 64._dp, 48._dp, 60._dp,&
2470 &                        0._dp, 32._dp,  2._dp, 36._dp,  0._dp,  0._dp,108._dp, 48._dp,  4._dp,120._dp,&
2471 &                        0._dp,  0._dp, 40._dp, 30._dp,  0._dp,  0._dp,  0._dp, 60._dp,120._dp,100._dp/),(/10,10/))
2472    bpp(1,:,:)=b2pp(:,:)/1225._dp
2473    write(std_out,*) "warning: this test is only valid if f4of2_sla=0"
2474  endif
2475  if(lcor==3) then
2476 !   app(0,:,:)=one
2477 !   a2pp(:,:)=RESHAPE((/ 100.0/1225.0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2478 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2479 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2480 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2481 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2482 !&                        0 ,0, 0, 0, 0, 100.0/1225.0, 0, 0, 0, 0, 0, 0, 0, 0,&
2483 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2484 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2485 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2486 !&                        0 ,0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,&
2487 !&                        0 ,0, 1,-1,-1, 1, 0, 0, 0, 0, 0, 0, 0, 0,&
2488 !&                        0 ,0,-1, 1, 1,-1, 0, 0, 0, 0, 0, 0, 0, 0,&
2489 !&                        0 ,0,-1, 1, 1,-1, 0, 0, 0, 0, 0, 0, 0, 0,&
2490 !&                        0 ,0, 1,-1,-1, 1, 0, 0, 0, 0, 0, 0, 0, 0,&/),(/14,14/))
2491    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,&
2492 &                        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,&
2493 &                        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,&
2494 &                        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,&
2495 &                        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,&
2496 &                        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,&
2497 &                        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,&
2498 &                        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,&
2499 &                        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,&
2500 &                        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,&
2501 &                        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,&
2502 &                        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,&
2503 &                        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,&
2504 &                        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/))
2505    bpp(0,:,:)=b2pp(:,:)
2506    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,&
2507 &                      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,&
2508 &                       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,&
2509 &                        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,&
2510 &                        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,&
2511 &                        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,&
2512 &                        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,&
2513 &                        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,&
2514 &                        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,&
2515 &                        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,&
2516 &                        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,&
2517 &                        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,&
2518 &                        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,&
2519 &                        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/))
2520    bpp(1,:,:)=b2pp(:,:)/1225._dp
2521    write(std_out,*) "warning: only 5/2 5/2 elements are given and only for exchange and F2 !"
2522 !   app(1,:,:)=a2pp(:,:)
2523 !   app(1,:,:)=app(1,:,:)/25
2524 !
2525 !   bpp(0,:,:)=one
2526 !   b2pp(:,:)=RESHAPE((/0,0,1,2,3,4,&
2527 !&                      0,0,4,3,2,1,&
2528 !&                      1,4,1,2,2,0,&
2529 !&                      2,3,2,1,0,2,&
2530 !&                      3,2,2,0,1,2,&
2531 !&                      4,1,0,2,2,1 /),(/6,6/))
2532 !   bpp(1,:,:)=b2pp(:,:)
2533 !   bpp(1,:,:)=bpp(1,:,:)/25
2534  endif
2535 
2536  do kk=0,lcor
2537    do m1=1,tndim
2538      do m2=1,tndim
2539        !write(6,*) kk,m1,m2
2540        !write(6,*) "--",fk(kk),aklmlmp(kk,m1,m2)
2541        !write(6,*) "--",fk(kk),bklmlmp(kk,m1,m2)
2542        udens(m1,m2)=udens(m1,m2)+fk(kk)*app(kk,m1,m2)
2543        jdens(m1,m2)=jdens(m1,m2)+fk(kk)*bpp(kk,m1,m2)
2544      enddo
2545    enddo
2546  enddo
2547  write(message,'(2x,a,3x,14f10.4)') " Direct Interaction Matrix from Inglis tables (in the JMJ basis) "
2548  call wrtout(std_out,  message,'COLL')
2549  write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,tndim,1)
2550  call wrtout(std_out,  message,'COLL')
2551  do m1=1,tndim
2552    write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2),m2=1,tndim,1)
2553    call wrtout(std_out,  message,'COLL')
2554  enddo
2555 
2556  write(message,'(a,2x,a,3x,14f10.4)') ch10," Exchange Interaction Matrix from Inglis tables (in the JMJ basis) "
2557  call wrtout(std_out,  message,'COLL')
2558  write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,tndim,1)
2559  call wrtout(std_out,  message,'COLL')
2560  do m1=1,tndim
2561    write(message,'(2x,i4,3x,14f10.4)') m1,(jdens(m1,m2),m2=1,tndim,1)
2562    call wrtout(std_out,  message,'COLL')
2563  enddo
2564 
2565  write(message,'(a,2x,a,3x,14f10.4)') ch10, " Density Density interactions from Inglis tables (in the JMJ basis) "
2566  call wrtout(std_out,  message,'COLL')
2567  write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,tndim,1)
2568  call wrtout(std_out,  message,'COLL')
2569  do m1=1,tndim
2570    write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2)-jdens(m1,m2),m2=1,tndim,1)
2571    call wrtout(std_out,  message,'COLL')
2572  enddo
2573 
2574 
2575  ABI_DEALLOCATE(jdens)
2576  ABI_DEALLOCATE(udens)
2577  ABI_DEALLOCATE(app)
2578  ABI_DEALLOCATE(bpp)
2579  ABI_DEALLOCATE(a2pp)
2580  ABI_DEALLOCATE(b2pp)
2581 
2582  end subroutine udens_inglis_hu

ABINIT/udens_slatercondon_hu [ Functions ]

[ Top ] [ 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-2018 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

2081 subroutine udens_slatercondon_hu(fk,lcor)
2082 
2083  use defs_basis
2084  use m_errors
2085  use m_abicore
2086 
2087 !This section has been created automatically by the script Abilint (TD).
2088 !Do not modify the following lines by hand.
2089 #undef ABI_FUNC
2090 #define ABI_FUNC 'udens_slatercondon_hu'
2091 !End of the abilint section
2092 
2093  implicit none
2094 
2095 !Arguments ---------------------------------------------
2096 !scalars
2097  integer,intent(in) :: lcor
2098  real(dp), intent(in) :: fk(0:lcor)
2099 
2100 !Local variables ---------------------------------------
2101 !scalars
2102  character(len=500) :: message
2103  integer :: m1,m2,kk
2104 !arrays
2105  real(dp), allocatable :: aklmlmp(:,:,:)
2106  real(dp), allocatable :: bklmlmp(:,:,:)
2107  real(dp), allocatable :: udens(:,:)
2108  real(dp), allocatable :: jdens(:,:)
2109 
2110 !*********************************************************************
2111  ABI_ALLOCATE(aklmlmp,(0:lcor,-lcor:lcor,-lcor:lcor)) ! k,m,m'
2112  ABI_ALLOCATE(bklmlmp,(0:lcor,-lcor:lcor,-lcor:lcor)) ! k,m,m'
2113  ABI_ALLOCATE(udens,(-lcor:lcor,-lcor:lcor)) ! m,m'
2114  ABI_ALLOCATE(jdens,(-lcor:lcor,-lcor:lcor)) ! m,m'
2115 ! k=2*(lcor)
2116  aklmlmp=zero
2117  bklmlmp=zero
2118  udens=zero
2119  jdens=zero
2120  if(lcor==0) then
2121    aklmlmp(0,0,0)=1
2122 !
2123    bklmlmp(0,0,0)=0
2124  else if(lcor==1) then
2125    aklmlmp(0, :, :)=1
2126    aklmlmp(1, 1, 1)= one/25._dp
2127    aklmlmp(1,-1,-1)= one/25._dp
2128    aklmlmp(1,-1, 1)= one/25._dp
2129    aklmlmp(1, 1,-1)= one/25._dp
2130    aklmlmp(1, 1, 0)=-two/25._dp
2131    aklmlmp(1,-1, 0)=-two/25._dp
2132    aklmlmp(1, 0,-1)=-two/25._dp
2133    aklmlmp(1, 0, 1)=-two/25._dp
2134    aklmlmp(1, 0, 0)= four/25._dp
2135 !
2136    bklmlmp(0, 1, 1)= one
2137    bklmlmp(0,-1,-1)= one
2138    bklmlmp(0, 0, 0)= one
2139    bklmlmp(1, 1, 1)= one/25._dp
2140    bklmlmp(1,-1,-1)= one/25._dp
2141    bklmlmp(1,-1,+1)= six/25._dp
2142    bklmlmp(1,+1,-1)= six/25._dp
2143    bklmlmp(1,+1, 0)= three/25._dp
2144    bklmlmp(1,-1, 0)= three/25._dp
2145    bklmlmp(1, 0,-1)= three/25._dp
2146    bklmlmp(1, 0, 1)= three/25._dp
2147    bklmlmp(1, 0, 0)= four/25._dp
2148  else if(lcor==2) then
2149    aklmlmp(0, :, :)=1
2150    aklmlmp(1, 2, 2)=  four / 49._dp
2151    aklmlmp(1, 2, 1)=  -two / 49._dp
2152    aklmlmp(1, 2, 0)= -four / 49._dp
2153    aklmlmp(1, 2,-1)=  -two / 49._dp
2154    aklmlmp(1, 2,-2)=  four / 49._dp
2155    aklmlmp(1, 1, 2)=  -two / 49._dp
2156    aklmlmp(1, 1, 1)=   one / 49._dp
2157    aklmlmp(1, 1, 0)=   two / 49._dp
2158    aklmlmp(1, 1,-1)=   one / 49._dp
2159    aklmlmp(1, 1,-2)=  -two / 49._dp
2160    aklmlmp(1, 0, 2)= -four / 49._dp
2161    aklmlmp(1, 0, 1)=   two / 49._dp
2162    aklmlmp(1, 0, 0)=  four / 49._dp
2163    aklmlmp(1, 0,-1)=   two / 49._dp
2164    aklmlmp(1, 0,-2)= -four / 49._dp
2165    aklmlmp(1,-1, 2)=  -two / 49._dp
2166    aklmlmp(1,-1, 1)=   one / 49._dp
2167    aklmlmp(1,-1, 0)=   two / 49._dp
2168    aklmlmp(1,-1,-1)=   one / 49._dp
2169    aklmlmp(1,-1,-2)=  -two / 49._dp
2170    aklmlmp(1,-2, 2)=  four / 49._dp
2171    aklmlmp(1,-2, 1)=  -two / 49._dp
2172    aklmlmp(1,-2, 0)= -four / 49._dp
2173    aklmlmp(1,-2,-1)=  -two / 49._dp
2174    aklmlmp(1,-2,-2)=  four / 49._dp
2175 
2176    aklmlmp(2, 2, 2)=    one  / 441._dp
2177    aklmlmp(2, 2, 1)=  -four  / 441._dp
2178    aklmlmp(2, 2, 0)=    six  / 441._dp
2179    aklmlmp(2, 2,-1)=  -four  / 441._dp
2180    aklmlmp(2, 2,-2)=    one  / 441._dp
2181    aklmlmp(2, 1, 2)=  -four  / 441._dp
2182    aklmlmp(2, 1, 1)=  16._dp / 441._dp
2183    aklmlmp(2, 1, 0)= -24._dp / 441._dp
2184    aklmlmp(2, 1,-1)=  16._dp / 441._dp
2185    aklmlmp(2, 1,-2)=  -four  / 441._dp
2186    aklmlmp(2, 0, 2)=    six  / 441._dp
2187    aklmlmp(2, 0, 1)= -24._dp / 441._dp
2188    aklmlmp(2, 0, 0)=  36._dp / 441._dp
2189    aklmlmp(2, 0,-1)= -24._dp / 441._dp
2190    aklmlmp(2, 0,-2)=    six  / 441._dp
2191    aklmlmp(2,-1, 2)=  -four  / 441._dp
2192    aklmlmp(2,-1, 1)=  16._dp / 441._dp
2193    aklmlmp(2,-1, 0)= -24._dp / 441._dp
2194    aklmlmp(2,-1,-1)=  16._dp / 441._dp
2195    aklmlmp(2,-1,-2)=  -four  / 441._dp
2196    aklmlmp(2,-2, 2)=    one  / 441._dp
2197    aklmlmp(2,-2, 1)=  -four  / 441._dp
2198    aklmlmp(2,-2, 0)=    six  / 441._dp
2199    aklmlmp(2,-2,-1)=  -four  / 441._dp
2200    aklmlmp(2,-2,-2)=    one  / 441._dp
2201    !do m1=lcor,-lcor,-1
2202    ! do m2=lcor,-lcor,-1
2203    !   write(6,*) m1,m2,aklmlmp(2,m1,m2)
2204    ! enddo
2205    !enddo
2206 
2207    !write(message,'(2x,a,3x,14f10.4)') " Slater aklmlmp(2,:,:)"
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,(aklmlmp(2,m1,m2),m2=-lcor,lcor,1)
2213    !  call wrtout(std_out,  message,'COLL')
2214    !enddo
2215 
2216    bklmlmp(0, 2, 2)=1
2217    bklmlmp(0, 1, 1)=1
2218    bklmlmp(0, 0, 0)=1
2219    bklmlmp(0,-1,-1)=1
2220    bklmlmp(0,-2,-2)=1
2221    bklmlmp(1, 2, 2)= four / 49._dp
2222    bklmlmp(1, 2, 1)=  six / 49._dp
2223    bklmlmp(1, 2, 0)= four / 49._dp
2224    bklmlmp(1, 2,-1)=  zero
2225    bklmlmp(1, 2,-2)=  zero
2226    bklmlmp(1, 1, 2)=  six / 49._dp
2227    bklmlmp(1, 1, 1)=  one / 49._dp
2228    bklmlmp(1, 1, 0)=  one / 49._dp
2229    bklmlmp(1, 1,-1)=  six / 49._dp
2230    bklmlmp(1, 1,-2)=  zero
2231    bklmlmp(1, 0, 2)= four / 49._dp
2232    bklmlmp(1, 0, 1)=  one / 49._dp
2233    bklmlmp(1, 0, 0)= four / 49._dp
2234    bklmlmp(1, 0,-1)=  one / 49._dp
2235    bklmlmp(1, 0,-2)= four / 49._dp
2236    bklmlmp(1,-1, 2)=  zero
2237    bklmlmp(1,-1, 1)=  six / 49._dp
2238    bklmlmp(1,-1, 0)=  one / 49._dp
2239    bklmlmp(1,-1,-1)=  one / 49._dp
2240    bklmlmp(1,-1,-2)=  six / 49._dp
2241    bklmlmp(1,-2, 2)=  zero
2242    bklmlmp(1,-2, 1)=  zero
2243    bklmlmp(1,-2, 0)= four / 49._dp
2244    bklmlmp(1,-2,-1)=  six / 49._dp
2245    bklmlmp(1,-2,-2)= four / 49._dp
2246 
2247    bklmlmp(2, 2, 2)=   one  / 441._dp
2248    bklmlmp(2, 2, 1)=  five  / 441._dp
2249    bklmlmp(2, 2, 0)= 15._dp / 441._dp
2250    bklmlmp(2, 2,-1)= 35._dp / 441._dp
2251    bklmlmp(2, 2,-2)= 70._dp / 441._dp
2252    bklmlmp(2, 1, 2)=  five  / 441._dp
2253    bklmlmp(2, 1, 1)= 16._dp / 441._dp
2254    bklmlmp(2, 1, 0)= 30._dp / 441._dp
2255    bklmlmp(2, 1,-1)= 40._dp / 441._dp
2256    bklmlmp(2, 1,-2)= 35._dp / 441._dp
2257    bklmlmp(2, 0, 2)= 15._dp / 441._dp
2258    bklmlmp(2, 0, 1)= 30._dp / 441._dp
2259    bklmlmp(2, 0, 0)= 36._dp / 441._dp
2260    bklmlmp(2, 0,-1)= 30._dp / 441._dp
2261    bklmlmp(2, 0,-2)= 15._dp / 441._dp
2262    bklmlmp(2,-1, 2)= 35._dp / 441._dp
2263    bklmlmp(2,-1, 1)= 40._dp / 441._dp
2264    bklmlmp(2,-1, 0)= 30._dp / 441._dp
2265    bklmlmp(2,-1,-1)= 16._dp / 441._dp
2266    bklmlmp(2,-1,-2)=  five  / 441._dp
2267    bklmlmp(2,-2, 2)= 70._dp / 441._dp
2268    bklmlmp(2,-2, 1)= 35._dp / 441._dp
2269    bklmlmp(2,-2, 0)= 15._dp / 441._dp
2270    bklmlmp(2,-2,-1)=  five  / 441._dp
2271    bklmlmp(2,-2,-2)=   one  / 441._dp
2272 
2273    !write(message,'(2x,a,3x,14f10.4)') " Slater bklmlmp(2,:,:)"
2274    !call wrtout(std_out,  message,'COLL')
2275    !write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1)
2276    !call wrtout(std_out,  message,'COLL')
2277    !do m1=-lcor,lcor,1
2278    !  write(message,'(2x,i4,3x,14f10.4)') m1,(bklmlmp(2,m1,m2),m2=-lcor,lcor,1)
2279    !  call wrtout(std_out,  message,'COLL')
2280    !enddo
2281          ! if f4of2_sla: these data agree with the explicit calculation
2282  else if(lcor==3) then
2283  endif
2284 
2285  do kk=0,lcor
2286    do m1=-lcor,lcor,1
2287      do m2=-lcor,lcor,1
2288        !write(6,*) kk,m1,m2
2289        !write(6,*) "--",fk(kk),aklmlmp(kk,m1,m2)
2290        !write(6,*) "--",fk(kk),bklmlmp(kk,m1,m2)
2291        udens(m1,m2)=udens(m1,m2)+fk(kk)*aklmlmp(kk,m1,m2)
2292        jdens(m1,m2)=jdens(m1,m2)+fk(kk)*bklmlmp(kk,m1,m2)
2293      enddo
2294    enddo
2295  enddo
2296  write(message,'(2x,a,3x,14f10.4)') " Direct Interaction Matrix from Slater tables (in the Ylm basis) "
2297  call wrtout(std_out,  message,'COLL')
2298  write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1)
2299  call wrtout(std_out,  message,'COLL')
2300  do m1=-lcor,lcor,1
2301    write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2),m2=-lcor,lcor,1)
2302    call wrtout(std_out,  message,'COLL')
2303  enddo
2304 
2305  write(message,'(a,2x,a,3x,14f10.4)') ch10," Exchange Interaction Matrix from Slater tables (in the Ylm basis) "
2306  call wrtout(std_out,  message,'COLL')
2307  write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1)
2308  call wrtout(std_out,  message,'COLL')
2309  do m1=-lcor,lcor,1
2310    write(message,'(2x,i4,3x,14f10.4)') m1,(jdens(m1,m2),m2=-lcor,lcor,1)
2311    call wrtout(std_out,  message,'COLL')
2312  enddo
2313 
2314  write(message,'(a,2x,a,3x,14f10.4)') ch10," Density density Interaction Matrix from Slater tables (in the Ylm basis) "
2315  call wrtout(std_out,  message,'COLL')
2316  write(message,'(2x,4x,14(2x,i8))') (m1,m1=-lcor,lcor,1)
2317  call wrtout(std_out,  message,'COLL')
2318  do m1=-lcor,lcor,1
2319    write(message,'(2x,i4,3x,14f10.4)') m1,(udens(m1,m2)-jdens(m1,m2),m2=-lcor,lcor,1)
2320    call wrtout(std_out,  message,'COLL')
2321  enddo
2322 
2323  ABI_DEALLOCATE(jdens)
2324  ABI_DEALLOCATE(udens)
2325  ABI_DEALLOCATE(aklmlmp)
2326  ABI_DEALLOCATE(bklmlmp)
2327 
2328  end subroutine udens_slatercondon_hu

ABINIT/vee_ylm2jmj_hu [ Functions ]

[ Top ] [ 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-2018 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

1915 subroutine vee_ylm2jmj_hu(lcor,mat_inp_c,mat_out_c,option)
1916 
1917  use defs_basis
1918  use m_errors
1919  use m_abicore
1920 
1921 !This section has been created automatically by the script Abilint (TD).
1922 !Do not modify the following lines by hand.
1923 #undef ABI_FUNC
1924 #define ABI_FUNC 'vee_ylm2jmj_hu'
1925 !End of the abilint section
1926 
1927  implicit none
1928 
1929 !Arguments ---------------------------------------------
1930 !scalars
1931  integer,intent(in) :: lcor,option
1932 !arrays
1933  complex(dpc),intent(in) :: mat_inp_c(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1))
1934  complex(dpc),intent(out) :: mat_out_c(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1))
1935 
1936 !Local variables ---------------------------------------
1937 !scalars
1938  integer :: ii,im,jc1,jj,jm,ll,ml1,ms1,gg,hh,gm,hm
1939  real(dp),parameter :: invsqrt2=one/sqrt2
1940  real(dp) :: invsqrt2lp1,xj,xmj
1941  complex(dpc) :: tmp2
1942  character(len=500) :: message
1943 !arrays
1944  integer, allocatable :: ind_msml(:,:)
1945  complex(dpc),allocatable :: mlms2jmj(:,:)
1946 
1947 !*********************************************************************
1948 
1949  if (option/=1.and.option/=2) then
1950    message=' option=/1 and =/2 !'
1951    MSG_BUG(message)
1952  end if
1953 
1954  if(option==1) then
1955    write(message,'(3a)') ch10,&
1956 &   "matrix in |l,s,m_l,m_s> basis is changed into |l,s,j,m_j> basis"
1957    call wrtout(std_out,message)
1958  else if(option==2) then
1959    write(message,'(3a)') ch10,&
1960 &   "matrix in |l,s,j,m_j> basis is changed into |l,s,m_l,m_s> basis"
1961    call wrtout(std_out,message)
1962  end if
1963 
1964 !--------------- Built indices + allocations
1965  ll=lcor
1966  ABI_ALLOCATE(mlms2jmj,(2*(2*ll+1),2*(2*ll+1)))
1967  mlms2jmj=czero
1968  ABI_ALLOCATE(ind_msml,(2,-ll:ll))
1969  mlms2jmj=czero
1970  jc1=0
1971  do ms1=1,2
1972    do ml1=-ll,ll
1973      jc1=jc1+1
1974      ind_msml(ms1,ml1)=jc1
1975    end do
1976  end do
1977 
1978 !--------------- built mlms2jmj
1979 !do jj=ll,ll+1    ! the physical value of j are ll-0.5,ll+0.5
1980 !xj(jj)=jj-0.5
1981  if(ll==0)then
1982    message=' ll should not be equal to zero !'
1983    MSG_BUG(message)
1984  end if
1985  jc1=0
1986  invsqrt2lp1=one/sqrt(float(2*lcor+1))
1987  do jj=ll,ll+1
1988    xj=float(jj)-half !  xj is in {ll-0.5, ll+0.5}
1989    do jm=-jj,jj-1
1990      xmj=float(jm)+half  ! xmj is in {-xj,xj}
1991      jc1=jc1+1           ! Global index for JMJ
1992      if(nint(xj+0.5)==ll+1) then  ! if xj=ll+0.5
1993        if(nint(xmj+0.5)==ll+1)  then
1994          mlms2jmj(ind_msml(2,ll),jc1)=1.0   !  J=L+0.5 and m_J=L+0.5
1995        else if(nint(xmj-0.5)==-ll-1) then
1996          mlms2jmj(ind_msml(1,-ll),jc1)=1.0   !  J=L+0.5 and m_J=-L-0.5
1997        else
1998          mlms2jmj(ind_msml(2,nint(xmj-0.5)),jc1)=invsqrt2lp1*(sqrt(float(ll)+xmj+0.5))
1999          mlms2jmj(ind_msml(1,nint(xmj+0.5)),jc1)=invsqrt2lp1*(sqrt(float(ll)-xmj+0.5))
2000        end if
2001      end if
2002      if(nint(xj+0.5)==ll) then  ! if xj=ll-0.5
2003        mlms2jmj(ind_msml(1,nint(xmj+0.5)),jc1)=invsqrt2lp1*(sqrt(float(ll)+xmj+0.5))
2004        mlms2jmj(ind_msml(2,nint(xmj-0.5)),jc1)=-invsqrt2lp1*(sqrt(float(ll)-xmj+0.5))
2005      end if
2006    end do
2007  end do
2008  write(message,'(3a)') ch10,"Matrix to go from |M_L,M_S> to |J,M_J>"
2009  call wrtout(std_out,message,"COLL")
2010  do im=1,2*(ll*2+1)
2011    write(message,'(12(1x,18(1x,"(",f7.3,",",f7.3,")")))') (mlms2jmj(im,jm),jm=1,2*(ll*2+1))
2012    call wrtout(std_out,message,"COLL")
2013  end do
2014 
2015 !--------------- compute change of basis
2016  do jm=1,2*(2*ll+1)
2017    do im=1,2*(2*ll+1)
2018      do hm=1,2*(2*ll+1)
2019        do gm=1,2*(2*ll+1)
2020          tmp2=czero
2021          do gg=1,2*(2*ll+1)
2022            do hh=1,2*(2*ll+1)
2023              do ii=1,2*(2*ll+1)
2024                do jj=1,2*(2*ll+1)
2025                  if(option==1) then
2026                    tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*CONJG(mlms2jmj(ii,im))*(mlms2jmj(jj,jm))&
2027 &                                                  *CONJG(mlms2jmj(gg,gm))*(mlms2jmj(hh,hm))
2028                  else if(option==2) then
2029 !                   tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*CONJG(mlms2jmj(ii,im))*(mlms2jmj(jj,jm))& ! inv=t*
2030 !&                                                  *CONJG(mlms2jmj(gg,gm))*(mlms2jmj(hh,hm)) ! inv=t*
2031                    tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*(mlms2jmj(im,ii))*CONJG(mlms2jmj(jm,jj))& ! inv=t*
2032 &                                                  *(mlms2jmj(gm,gg))*CONJG(mlms2jmj(hm,hh)) ! inv=t*
2033                  end if
2034                end do
2035              end do
2036            end do
2037          end do
2038          mat_out_c(gm,im,hm,jm)=tmp2
2039        end do
2040      end do
2041    end do
2042  end do
2043  ABI_DEALLOCATE(mlms2jmj)
2044  ABI_DEALLOCATE(ind_msml)
2045 
2046  end subroutine vee_ylm2jmj_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

PARENTS

      m_dmft

CHILDREN

      wrtout

SOURCE

329 subroutine destroy_hu(hu,ntypat,t2g)
330 
331  use defs_basis
332  use m_crystal, only : crystal_t
333 
334 !This section has been created automatically by the script Abilint (TD).
335 !Do not modify the following lines by hand.
336 #undef ABI_FUNC
337 #define ABI_FUNC 'destroy_hu'
338 !End of the abilint section
339 
340  implicit none
341 
342 !Arguments ------------------------------------
343 !scalars
344  integer, intent(in) :: ntypat
345  type(hu_type),intent(inout) :: hu(ntypat)
346  integer, intent(in) :: t2g
347 
348 !Local variables-------------------------------
349  integer :: itypat
350 
351 ! *********************************************************************
352 
353  do itypat=1,ntypat
354 !  if ( allocated(hu(itypat)%vee) )  deallocate(hu(itypat)%vee)
355   if ( allocated(hu(itypat)%uqmc) )   then
356     ABI_DEALLOCATE(hu(itypat)%uqmc)
357   end if
358   if ( allocated(hu(itypat)%udens) )   then
359     ABI_DEALLOCATE(hu(itypat)%udens)
360   end if
361   if(t2g==1.and.hu(itypat)%lpawu==1) then
362     ABI_DEALLOCATE(hu(itypat)%vee)
363   else
364     hu(itypat)%vee => null()
365   endif
366  enddo
367 
368 end subroutine destroy_hu

m_hu/hu_type [ Types ]

[ Top ] [ m_hu ] [ Types ]

NAME

  hu_type

FUNCTION

  This structured datatype contains interaction matrices for the correlated subspace

SOURCE

68  type, public :: hu_type ! for each typat
69 
70   integer :: lpawu
71 
72   logical :: jmjbasis
73 
74   real(dp) :: upawu    ! => upaw
75 
76   real(dp) :: jpawu    ! => jpaw
77 
78   real(dp) :: f2_sla    ! => f2_sla
79 
80   real(dp) :: f4of2_sla    ! => f4of2_sla
81 
82   real(dp) :: f6of2_sla    ! => f6of2_sla
83 
84   logical :: jpawu_zero  ! true if all jpawu are zero
85                          ! false if one of the jpaw is not zero
86 
87   real(dp), pointer :: vee(:,:,:,:) => null() ! => vee
88 
89   real(dp), allocatable :: uqmc(:)
90 
91   real(dp), allocatable :: udens(:,:)
92 
93  end type hu_type
94 
95 !----------------------------------------------------------------------
96 
97 
98 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

PARENTS

      m_dmft

CHILDREN

      wrtout

SOURCE

123 subroutine init_hu(cryst_struc,pawtab,hu,t2g)
124 
125  use defs_basis
126  use m_crystal, only : crystal_t
127 
128 !This section has been created automatically by the script Abilint (TD).
129 !Do not modify the following lines by hand.
130 #undef ABI_FUNC
131 #define ABI_FUNC 'init_hu'
132 !End of the abilint section
133 
134  implicit none
135 
136 !Arguments ------------------------------------
137 !type
138  type(crystal_t),intent(in) :: cryst_struc
139  type(pawtab_type), target, intent(in)  :: pawtab(cryst_struc%ntypat)
140  type(hu_type), intent(inout) :: hu(cryst_struc%ntypat)
141  integer, intent(in) :: t2g
142 !Local variables ------------------------------------
143  integer :: itypat,i,ij,ij1,ij2,j,lpawu,ms,ms1,m,m1,ndim
144  integer :: ns,ns1,n,n1
145  integer, allocatable :: xij(:,:)
146  real(dp) :: xtemp
147  character(len=500) :: message
148 !************************************************************************
149  write(message,'(2a)') ch10,"  == Compute Interactions for DMFT"
150  call wrtout(std_out,message,'COLL')
151 
152  xtemp=zero
153 
154 ! ====================================
155 !  Compute hu(iatom)%uqmc from vee
156 ! ====================================
157  hu(1)%jpawu_zero=.true.
158  do itypat=1,cryst_struc%ntypat
159    hu(itypat)%lpawu=pawtab(itypat)%lpawu
160    hu(itypat)%jmjbasis=.false.
161    if(t2g==1.and.hu(itypat)%lpawu==2) hu(itypat)%lpawu=1
162    lpawu=hu(itypat)%lpawu
163    if(lpawu.ne.-1) then
164      hu(itypat)%upawu=pawtab(itypat)%upawu
165      hu(itypat)%jpawu=pawtab(itypat)%jpawu
166 
167    ! The if below are necessary for an unknown reason with the NAG  compiler.
168     if(pawtab(itypat)%f4of2_sla>=0) then
169       hu(itypat)%f4of2_sla = pawtab(itypat)%f4of2_sla
170     else if(pawtab(itypat)%f4of2_sla==0) then
171       hu(itypat)%f4of2_sla = 0_dp
172     else if(pawtab(itypat)%f4of2_sla<0) then
173       hu(itypat)%f4of2_sla = pawtab(itypat)%f4of2_sla
174     endif
175     if(pawtab(itypat)%f6of2_sla>=0) then
176       hu(itypat)%f6of2_sla = pawtab(itypat)%f6of2_sla
177     else if(pawtab(itypat)%f6of2_sla==0) then
178       hu(itypat)%f6of2_sla = 0_dp
179     else if(pawtab(itypat)%f6of2_sla<0) then
180       hu(itypat)%f6of2_sla = pawtab(itypat)%f6of2_sla
181     endif
182 
183 !    This is copied from pawpuxinit: it would be better not to duplicate
184 !    these lines.
185      if(lpawu==0) then
186        hu(itypat)%f2_sla=zero
187      else if(lpawu==1) then
188        hu(itypat)%f2_sla=pawtab(itypat)%jpawu*5._dp
189      else if(lpawu==2) then
190        hu(itypat)%f2_sla=pawtab(itypat)%jpawu*14._dp/(One+hu(itypat)%f4of2_sla)
191      else if(lpawu==3) then
192        hu(itypat)%f2_sla=pawtab(itypat)%jpawu*6435._dp/(286._dp+&
193 &       195._dp*hu(itypat)%f4of2_sla+250._dp*hu(itypat)%f6of2_sla)
194      endif
195 
196      if(hu(itypat)%jpawu>tol4) hu(1)%jpawu_zero=.false.
197      ndim=2*lpawu+1
198 !     ndim1=2*hu(itypat)%lpawu+1
199      write(message,'(2a,i4)')  ch10,'  -------> For Correlated Species', itypat
200      call wrtout(std_out,  message,'COLL')
201 !     allocate(hu(itypat)%vee(ndim,ndim,ndim,ndim))
202      ABI_ALLOCATE(hu(itypat)%uqmc,(ndim*(2*ndim-1)))
203      ABI_ALLOCATE(hu(itypat)%udens,(2*ndim,2*ndim))
204      ABI_ALLOCATE(xij,(2*ndim,2*ndim))
205      if(t2g==0) then
206        hu(itypat)%vee => pawtab(itypat)%vee
207 !   t2g case begin
208      else if(t2g==1.and.hu(itypat)%lpawu==1) then
209        ABI_ALLOCATE(hu(itypat)%vee,(ndim,ndim,ndim,ndim))
210        n=0
211        do m=1,5
212          if((m/=1.and.m/=2.and.m/=4)) cycle
213          n=n+1
214          ns=0
215          do ms=1,5
216            if((ms/=1.and.ms/=2.and.ms/=4)) cycle
217            ns=ns+1
218            n1=0
219            do m1=1,5
220              if((m1/=1.and.m1/=2.and.m1/=4)) cycle
221              n1=n1+1
222              ns1=0
223              do ms1=1,5
224                if((ms1/=1.and.ms1/=2.and.ms1/=4)) cycle
225                ns1=ns1+1
226                hu(itypat)%vee(n,ns,n1,ns1)=pawtab(itypat)%vee(m,ms,m1,ms1)
227              enddo
228            enddo
229          enddo
230        enddo
231 !   t2g case end
232      endif
233 
234      hu(itypat)%udens=zero
235      ij=0
236      do ms=1,2*ndim-1
237          xij(ms,ms)=0
238        do ms1=ms+1,2*ndim
239          ij=ij+1
240          xij(ms,ms1)=ij
241          xij(ms1,ms)=ij
242          if(ms<=ndim.and.ms1>ndim) then
243            m1 = ms1 - ndim
244            m  = ms
245            hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)
246            hu(itypat)%udens(ms,ms1)= hu(itypat)%vee(m,m1,m,m1)
247            hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1)
248          else if(ms<=ndim.and.ms1<=ndim) then
249            m1 = ms1
250            m  = ms
251            hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m)
252            hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij)
253            hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1)
254          else
255            m1 = ms1 - ndim
256            m  = ms  - ndim
257            hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m)
258            hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij)
259            hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1)
260          endif
261        enddo
262      enddo
263      xij(2*ndim,2*ndim)=0
264      write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation "
265      call wrtout(std_out,  message,'COLL')
266      write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim)
267      call wrtout(std_out,  message,'COLL')
268 !     xtemp1b=0.d0
269 ! ====================================
270 !  Print hu(iatom)%uqmc
271 ! ====================================
272      ij1=-10
273      ij2=-10
274      ij=0
275      do i=1,2*ndim
276        do j=i+1,2*ndim
277          ij=ij+1
278          if(j==i+1) ij1=ij
279          if(j==2*ndim) ij2=ij
280        enddo
281 !       write(std_out,*) itypat
282 !       do m=1,i
283 !        write(std_out,*) i,m
284 !        write(std_out,*) xij(i,m)
285 !        write(std_out,*) ij1,ij2
286 !       enddo
287        if(i==1)               write(message,'(i3,14f7.3)') &
288 &                              i,xtemp, (hu(itypat)%uqmc(m),m=ij1,ij2)
289        if(i/=2*ndim.and.i/=1) write(message,'(i3,14f7.3)') i, &
290 &        (hu(itypat)%uqmc(xij(i,m)), m=1,i-1),xtemp, (hu(itypat)%uqmc(m),m=ij1,ij2)
291        if(i==2*ndim)          write(message,'(i3,14f7.3)') i, &
292 &                  (hu(itypat)%uqmc(xij(i,m)), m=1,i-1),xtemp
293        call wrtout(std_out,  message,'COLL')
294      enddo
295        write(message,'(5x,a)') "--------------------------------------------------------"
296        call wrtout(std_out,  message,'COLL')
297      ABI_DEALLOCATE(xij)
298    else
299      hu(itypat)%upawu=zero
300      hu(itypat)%jpawu=zero
301 !     allocate(hu(itypat)%vee(0,0,0,0))
302    endif
303  enddo ! itypat
304 
305 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

393 subroutine print_hu(hu,ntypat,prtopt)
394 
395  use defs_basis
396  use m_crystal, only : crystal_t
397 
398 !This section has been created automatically by the script Abilint (TD).
399 !Do not modify the following lines by hand.
400 #undef ABI_FUNC
401 #define ABI_FUNC 'print_hu'
402 !End of the abilint section
403 
404  implicit none
405 
406 !Arguments ------------------------------------
407 !type
408  integer, intent(in):: ntypat
409  type(hu_type),intent(in) :: hu(ntypat)
410  integer :: prtopt
411 
412 !Local variables-------------------------------
413  integer :: itypat
414  integer :: lpawu,ms,ms1,m,ndim
415  character(len=500) :: message
416 ! *********************************************************************
417 
418  do itypat = 1 , ntypat
419    lpawu=hu(itypat)%lpawu
420    if(lpawu/=-1) then
421      ndim=2*lpawu+1
422      write(message,'(2a,i4)')  ch10,'  -------> For Correlated species'
423      call wrtout(std_out,  message,'COLL')
424      if(prtopt==0) then
425        write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in Slm basis "
426      else if(prtopt==1) then
427        write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in diagonal basis"
428      else if(prtopt==2) then
429        write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in Ylm basis"
430      else if(prtopt==3) then
431        write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation in JMJ basis"
432      endif
433      call wrtout(std_out,  message,'COLL')
434      write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim)
435      call wrtout(std_out,  message,'COLL')
436        do ms=1,2*ndim
437           write(message,'(i3,14f7.3)') &
438 &          ms, (hu(itypat)%udens(ms,ms1),ms1=1,2*ndim)
439           call wrtout(std_out,  message,'COLL')
440        enddo
441        write(message,'(5x,a)') "--------------------------------------------------------"
442        call wrtout(std_out,  message,'COLL')
443    endif ! lpawu/=1
444  enddo ! ntypat
445 
446 
447 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

1057 subroutine printvee_hu(ndim,vee,prtopt,basis,upawu,f2)
1058 
1059  use defs_basis
1060  use m_crystal, only : crystal_t
1061 
1062 !This section has been created automatically by the script Abilint (TD).
1063 !Do not modify the following lines by hand.
1064 #undef ABI_FUNC
1065 #define ABI_FUNC 'printvee_hu'
1066 !End of the abilint section
1067 
1068  implicit none
1069 
1070 !Arguments ------------------------------------
1071 !type
1072  integer,intent(in) :: ndim,prtopt
1073  real(dp), intent(in) :: vee(ndim,ndim,ndim,ndim)
1074  real(dp), optional, intent(in) :: upawu,f2
1075  character(len=*), intent(in) :: basis
1076 
1077 !Local variables-------------------------------
1078  integer :: abcomp,m1,m2,mi,mj,mk,ml
1079  real(dp),allocatable :: b0(:,:)
1080  real(dp),allocatable :: a2pp(:,:)
1081  real(dp),allocatable :: b2pp(:,:)
1082  real(dp):: aver52,aver72,avernondiag,averall,averu,averj,averurestricted
1083  character(len=2000) :: message
1084 ! *********************************************************************
1085   write(message,'(5a)') ch10,&
1086 &  '  Coulomb interaction in the ', trim(basis),' basis'
1087   call wrtout(std_out,message,'COLL')
1088 
1089  if(prtopt>=2) then
1090 
1091    write(message,'(2a)')  ch10," <mi,mi|vee|mi mi> : U1"
1092    call wrtout(std_out,  message,'COLL')
1093    do mi=1,ndim
1094      write(message,'(4i4,3x,e10.3)')   mi,mi,mi,mi,vee(mi,mi,mi,mi)
1095      call wrtout(std_out,  message,'COLL')
1096    enddo
1097 
1098    write(message,'(2a)')  ch10," <mi,mj|vee|mi mj> : U2"
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,mj,mi,mj,vee(mi,mj,mi,mj)
1103        call wrtout(std_out,  message,'COLL')
1104      enddo
1105    enddo
1106 
1107    write(message,'(2a)')  ch10," <mi,mj|vee|mj mi> : J"
1108    call wrtout(std_out,  message,'COLL')
1109    do mi=1,ndim
1110      do mj=mi+1,ndim
1111        write(message,'(4i4,3x,e10.3)')   mi,mj,mj,mi,vee(mi,mj,mj,mi)
1112        call wrtout(std_out,  message,'COLL')
1113      enddo
1114    enddo
1115 
1116    write(message,'(2a)')  ch10," <mi,mi|vee|mj mj> : J"
1117    call wrtout(std_out,  message,'COLL')
1118    do mi=1,ndim
1119      do mj=mi+1,ndim
1120        write(message,'(4i4,3x,e10.3)')   mi,mi,mj,mj,vee(mi,mi,mj,mj)
1121        call wrtout(std_out,  message,'COLL')
1122      enddo
1123    enddo
1124 
1125    write(message,'(2a)')  ch10," vee is non zero also for"
1126    call wrtout(std_out,  message,'COLL')
1127    do mi=1,ndim
1128      do mj=1,ndim
1129        do mk=1,ndim
1130          do ml=1,ndim
1131             if((.not.(mi==mk.and.mj==ml)).and.(.not.(mi==ml.and.mj==mk)).and.&
1132 &               .not.(mi==mj.and.mk==ml)) then
1133               if(vee(mi,mj,mk,ml)>tol8) then
1134                 write(message,'(4i4,3x,e10.3)')   mi,mj,mk,ml,vee(mi,mj,mk,ml)
1135                 call wrtout(std_out,  message,'COLL')
1136               endif
1137             endif
1138          enddo
1139        enddo
1140      enddo
1141    enddo
1142    write(message,'(a)')  ch10
1143    call wrtout(std_out,  message,'COLL')
1144 
1145  endif
1146 
1147  if(prtopt>=1) then
1148 
1149    write(message,'(2x,a,3x,14f10.4)') "Um1m2=Vee(m1,m2,m1,m2)"
1150    call wrtout(std_out,  message,'COLL')
1151    write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim)
1152    call wrtout(std_out,  message,'COLL')
1153    do m1=1,ndim
1154      write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m1,m2),m2=1,ndim)
1155      call wrtout(std_out,  message,'COLL')
1156    enddo
1157    write(message,'(a)') ch10;  call wrtout(std_out,  message,'COLL')
1158 
1159    write(message,'(2x,a,3x,14f10.4)') "Jm1m2=Vee(m1,m2,m2,m1)"
1160    call wrtout(std_out,  message,'COLL')
1161    write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim)
1162    call wrtout(std_out,  message,'COLL')
1163    do m1=1,ndim
1164      write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m2,m1),m2=1,ndim)
1165      call wrtout(std_out,  message,'COLL')
1166    enddo
1167    write(message,'(a)') ch10;  call wrtout(std_out,  message,'COLL')
1168 
1169    write(message,'(2x,a,3x,14f10.4)') "Udens(m1,m2)"
1170    call wrtout(std_out,  message,'COLL')
1171    write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim)
1172    call wrtout(std_out,  message,'COLL')
1173    do m1=1,ndim
1174      write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1),m2=1,ndim)
1175      call wrtout(std_out,  message,'COLL')
1176    enddo
1177    write(message,'(a)') ch10;  call wrtout(std_out,  message,'COLL')
1178 
1179    if(prtopt>=3) then
1180      if(ndim==7) then
1181        averu=zero
1182        do m1=1,7
1183          do m2=1,7
1184            averu=vee(m1,m2,m1,m2)+averu
1185          enddo
1186        enddo
1187        averu=averu/49.d0
1188        averurestricted=zero
1189        do m1=1,7
1190          do m2=1,7
1191            if(m1/=m2) averurestricted=vee(m1,m2,m1,m2)+averurestricted
1192          enddo
1193        enddo
1194        averurestricted=averurestricted/42.d0
1195        averj=zero
1196        do m1=1,7
1197          do m2=1,7
1198            if(m1/=m2) averj=vee(m1,m2,m2,m1)+averj
1199          enddo
1200        enddo
1201        averj=averj/42
1202        averall=zero
1203        do m1=1,7
1204          do m2=1,7
1205            if(m1/=m2) averall=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+averall
1206          enddo
1207        enddo
1208        averall=averall/42
1209          !write(6,*) "averages Ylm U, U restricted,J, U-J",averu,averurestricted,averj,averall
1210 
1211      endif
1212 
1213      if(ndim==14.and.trim(basis)=='jmj') then
1214        aver52=zero
1215        do m1=1,6
1216          do m2=1,6
1217            aver52=vee(m1,m2,m1,m2)+aver52
1218          enddo
1219        enddo
1220        aver52=aver52/36.d0
1221        aver72=zero
1222        do m1=7,14
1223          do m2=7,14
1224             aver72=vee(m1,m2,m1,m2)+aver72
1225          enddo
1226        enddo
1227        aver72=aver72/64.d0
1228        avernondiag=zero
1229        do m1=1,6
1230          do m2=7,14
1231            avernondiag=vee(m1,m2,m1,m2)+avernondiag
1232          enddo
1233        enddo
1234        avernondiag=avernondiag/48.d0
1235        averall=zero
1236        do m1=1,14
1237          do m2=1,14
1238            averall=vee(m1,m2,m1,m2)+averall
1239          enddo
1240        enddo
1241        averall=averall/196.d0
1242          !write(6,*) "U averages",aver52,aver72,avernondiag,averall
1243 
1244 
1245 
1246        aver52=zero
1247        do m1=1,6
1248          do m2=1,6
1249            if(m1/=m2) aver52=vee(m1,m2,m2,m1)+aver52
1250          enddo
1251        enddo
1252        aver52=aver52/30.d0
1253        aver72=zero
1254        do m1=7,14
1255          do m2=7,14
1256            if(m1/=m2) aver72=vee(m1,m2,m2,m1)+aver72
1257          enddo
1258        enddo
1259        aver72=aver72/56.d0
1260        avernondiag=zero
1261        do m1=1,6
1262          do m2=7,14
1263            avernondiag=vee(m1,m2,m2,m1)+avernondiag
1264          enddo
1265        enddo
1266        avernondiag=avernondiag/48.d0
1267        averall=zero
1268        do m1=1,14
1269          do m2=1,14
1270            if(m1/=m2) averall=vee(m1,m2,m2,m1)+averall
1271          enddo
1272        enddo
1273        averall=averall/182.d0
1274          !write(6,*) "J averages",aver52,aver72,avernondiag,averall
1275 
1276 
1277 
1278 
1279        aver52=zero
1280        do m1=1,6
1281          do m2=1,6
1282            if(m1/=m2) aver52=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+aver52
1283          enddo
1284        enddo
1285        aver52=aver52/30.d0
1286        aver72=zero
1287        do m1=7,14
1288          do m2=7,14
1289            if(m1/=m2) aver72=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+aver72
1290          enddo
1291        enddo
1292        aver72=aver72/56.d0
1293        avernondiag=zero
1294        do m1=1,6
1295          do m2=7,14
1296            avernondiag=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+avernondiag
1297          enddo
1298        enddo
1299        avernondiag=avernondiag/48.d0
1300        averall=zero
1301        do m1=1,14
1302          do m2=1,14
1303            if(m1/=m2) averall=vee(m1,m2,m1,m2)-vee(m1,m2,m2,m1)+averall
1304          enddo
1305        enddo
1306        averall=averall/182.d0
1307        !write(6,*) "U-J averages",aver52,aver72,avernondiag,averall
1308 
1309      endif
1310    endif
1311 
1312    if (present(upawu)) then
1313      ABI_ALLOCATE(a2pp,(ndim,ndim))
1314      ABI_ALLOCATE(b2pp,(ndim,ndim))
1315      ABI_ALLOCATE(b0,(ndim,ndim))
1316 
1317 !     write(message,'(2x,a,3x,14f10.4)') "For check with respect to Slater's paper"
1318 !     call wrtout(std_out,  message,'COLL')
1319 !     ABI_ALLOCATE(f0,(ndim,ndim,ndim,ndim))
1320 !     write(message,'(2x,a,3x,14f10.4)') "Vee(m1,m2,m1,m2)-F0*ao(m1,m2)"
1321 !     call wrtout(std_out,  message,'COLL')
1322 !     write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim)
1323 !     call wrtout(std_out,  message,'COLL')
1324 !
1325 !     do m1=1,ndim
1326 !       write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m1,m2)-upawu,m2=1,ndim)
1327 !       call wrtout(std_out,  message,'COLL')
1328 !     enddo
1329 !
1330 !     f0=zero
1331 !     do m1=1,ndim
1332 !     f0(m1,m1,m1,m1)=upawu
1333 !     enddo
1334 !     write(message,'(a)') ch10
1335 !     call wrtout(std_out,  message,'COLL')
1336 !     write(message,'(2x,a,3x,14f10.4)') "Vee(m1,m2,m2,m1)-F0*b0(m1,m2)"
1337 !     call wrtout(std_out,  message,'COLL')
1338 !     write(message,'(2x,4x,14(2x,i8))') (m1,m1=1,ndim)
1339 !     call wrtout(std_out,  message,'COLL')
1340 !     do m1=1,ndim
1341 !       write(message,'(2x,i4,3x,14f10.4)') m1,(vee(m1,m2,m2,m1)-f0(m1,m2,m2,m1),m2=1,ndim)
1342 !       call wrtout(std_out,  message,'COLL')
1343 !     enddo
1344 !     ABI_DEALLOCATE(f0)
1345 
1346 
1347      b0=zero
1348      do m1=1,ndim
1349      b0(m1,m1)=upawu
1350      enddo
1351      abcomp=0
1352      if(ndim==3.and.present(f2).and.(trim(basis)=='slm')) then
1353        a2pp(:,:) = RESHAPE((/1,-2,1,-2,4,-2,-1,-2,-1/),(/3,3/))
1354        a2pp=a2pp/25*f2+upawu
1355        b2pp(:,:) = RESHAPE((/1,3,6,3,4,3,6,3,1/),(/3,3/))
1356        b2pp=b2pp/25*f2+b0
1357        abcomp=1
1358      else if(ndim==6.and.present(f2).and.(trim(basis)=='jmj')) then
1359        ABI_ALLOCATE(a2pp,(6,6))
1360        ABI_ALLOCATE(b2pp,(6,6))
1361        a2pp(:,:)=RESHAPE((/0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,-1,-1,1,0,&
1362 &       0,-1,1,1,-1,0,0,-1,1,1,-1,0,0,1,-1,-1,1/),(/6,6/))
1363        a2pp=a2pp/25*f2+upawu
1364        b2pp(:,:)=RESHAPE((/0,0,1,2,3,4,0,0,4,3,2,1,1,4,1,2,2,0,2,3,&
1365 &       2,1,0,2,3,2,2,0,1,2,4,1,0,2,2,1/),(/6,6/))
1366        b2pp=b2pp/25*f2+b0
1367        abcomp=1
1368      endif
1369      if(mod(ndim,3)==0.and.present(f2).and.abcomp==1) then
1370        write(message,'(2x,a)') "Exact result for Umm is"
1371        call wrtout(std_out,message,'COLL')
1372        do m1=1,ndim
1373          write(message,'(2x,i4,3x,14f10.4)') m1,(a2pp(m1,m2),m2=1,ndim)
1374          call wrtout(std_out,  message,'COLL')
1375        enddo
1376        write(message,'(a)') ch10;  call wrtout(std_out,  message,'COLL')
1377        write(message,'(2x,a,3x,14f10.4)') "Exact result for Jmm is"
1378        call wrtout(std_out,message,'COLL')
1379        do m1=1,ndim
1380          write(message,'(2x,i4,3x,14f10.4)') m1,(b2pp(m1,m2),m2=1,ndim)
1381          call wrtout(std_out,message,'COLL')
1382        enddo
1383        write(message,'(a)') ch10;  call wrtout(std_out,  message,'COLL')
1384      endif
1385      ABI_DEALLOCATE(a2pp)
1386      ABI_DEALLOCATE(b2pp)
1387      ABI_DEALLOCATE(b0)
1388    endif
1389 
1390  endif
1391 
1392 end subroutine printvee_hu

m_hu/reddd [ Functions ]

[ Top ] [ m_hu ] [ Functions ]

NAME

 reddd

FUNCTION

INPUTS

OUTPUT

PARENTS

CHILDREN

SOURCE

1517 function reddd(mi,ndim)
1518 
1519  use defs_basis
1520 
1521 !This section has been created automatically by the script Abilint (TD).
1522 !Do not modify the following lines by hand.
1523 #undef ABI_FUNC
1524 #define ABI_FUNC 'reddd'
1525 !End of the abilint section
1526 
1527  implicit none
1528 
1529 !Arguments ------------------------------------
1530 !scalars
1531  integer,intent(in) :: mi,ndim
1532  integer :: reddd
1533 ! *************************************************************************
1534 
1535  if(mi<ndim+1)  reddd=mi
1536  if(mi>=ndim+1) reddd=mi-ndim
1537 
1538 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.

PARENTS

      hubbard_one,qmc_prep_ctqmc

CHILDREN

      wrtout

SOURCE

 570 subroutine rotatevee_hu(cryst_struc,hu,nspinor,nsppol,pawprtvol,rot_mat,udens_atoms,rot_type)
 571 
 572  use defs_basis
 573  use m_errors
 574 
 575  use m_crystal,          only : crystal_t
 576 
 577 !This section has been created automatically by the script Abilint (TD).
 578 !Do not modify the following lines by hand.
 579 #undef ABI_FUNC
 580 #define ABI_FUNC 'rotatevee_hu'
 581 !End of the abilint section
 582 
 583  implicit none
 584 
 585 !Arguments ------------------------------------
 586 !type
 587  type(crystal_t),intent(in) :: cryst_struc
 588  integer, intent(in):: nsppol,nspinor,pawprtvol,rot_type
 589  type(hu_type),intent(in) :: hu(cryst_struc%ntypat)
 590  type(coeff2c_type),optional,intent(in) :: rot_mat(cryst_struc%natom,nsppol)
 591  type(coeff2_type),intent(inout) :: udens_atoms(cryst_struc%natom)
 592 
 593 !Local variables-------------------------------
 594  integer :: dim_vee,iatom,itypat
 595  integer :: lpawu,m1,m2,m3,m4,mi,mj,mk,ml,natom,ndim,tndim
 596  integer :: ms,ms1,nat_correl
 597  character(len=500) :: message
 598  character(len=30) :: basis_vee
 599  real(dp) :: xsum,xsum2,xsumnew,xsum2new,uaver
 600  complex(dpc),allocatable :: temp_mat(:,:)
 601  complex(dpc),allocatable :: temp_mat2(:,:)
 602  real(dp),allocatable :: veetemp(:,:,:,:),fk(:)
 603  complex(dpc),allocatable :: veeylm(:,:,:,:),veeslm(:,:,:,:)
 604  complex(dpc),allocatable :: veetmp(:,:,:),veetmp_slm(:,:,:),veetmp_ylm(:,:,:)
 605  complex(dpc),allocatable :: veejmj(:,:,:,:),veerotated(:,:,:,:),veenew(:,:,:,:)
 606  complex(dpc),allocatable :: veeylm2(:,:,:,:),veeslm2(:,:,:,:)
 607 ! *********************************************************************
 608  natom=cryst_struc%natom
 609 
 610 !================================================
 611 !  NSPINOR = 2 and J=0
 612 !================================================
 613  if(hu(1)%jpawu_zero.and.nspinor==2) then
 614 ! if(3==4) then
 615 
 616 !   call vee2udens_hu(hu,cryst_struc%ntypat,2)
 617    do iatom=1,natom
 618      itypat=cryst_struc%typat(iatom)
 619      lpawu=hu(itypat)%lpawu
 620      if(lpawu.ne.-1) then
 621        ndim=2*lpawu+1
 622        write(message,'(2a,i4)')  ch10,'  -------> For Correlated Species', itypat
 623        call wrtout(std_out,  message,'COLL')
 624     !   write(std_out,*)"ndim",ndim
 625     !   write(std_out,'(a,20e10.4)')"udens before vee2udensaomt", udens_atoms(1)%value
 626 !       write(std_out,*)"vee",hu(cryst_struc%typat(iatom))%vee
 627        call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,hu(cryst_struc%typat(iatom))%vee,"Slm")
 628     !   write(std_out,*)"udens after vee2udensatom", udens_atoms(1)%value
 629      endif
 630    enddo
 631   ! write(std_out,*)"udensafter after", udens_atoms(1)%value
 632 
 633 
 634 !================================================
 635 !  NSPINOR = 2 and J/=0
 636 !================================================
 637  else if (nspinor==2) then
 638    write(std_out,*)
 639    write(std_out,*)" == Rotate interaction in the level basis "
 640 
 641 
 642    !rot_type=0 ! keep original Slm basis
 643    !rot_type=2 ! rotation to the Ylm basis ! for tests
 644    !rot_type=4 ! use a rotation from Ylm basis to a new basis
 645    !rot_type=3 ! rotation to the JmJ Basis ! for tests
 646    !rot_type=1 ! use a rotation for diago of dmat,green, levels..
 647 
 648    do iatom=1,natom
 649      itypat=cryst_struc%typat(iatom)
 650      lpawu=hu(itypat)%lpawu
 651      if(lpawu.ne.-1) then
 652        ndim=2*lpawu+1
 653        if(pawprtvol>=3) then
 654 !         write(message,'(2a)')  ch10," VEE INPUT AVANT TRANSFORMATION"
 655 !         call wrtout(std_out,  message,'COLL')
 656 !         call printvee_hu(ndim,hu(itypat)%vee,1,'Slm')
 657        endif
 658        write(message,'(2a,i4)')  ch10,'  -------> For Correlated atom', iatom
 659        call wrtout(std_out,  message,'COLL')
 660        tndim=nspinor*ndim
 661        ABI_ALLOCATE(veejmj,(tndim,tndim,tndim,tndim))
 662        ABI_ALLOCATE(veeslm,(ndim,ndim,ndim,ndim))
 663        ABI_ALLOCATE(veetmp_slm,(ndim,ndim,4))
 664        ABI_ALLOCATE(veetmp_ylm,(ndim,ndim,4))
 665        ABI_ALLOCATE(veetmp,(ndim,ndim,4))
 666        ABI_ALLOCATE(veeslm2,(tndim,tndim,tndim,tndim))
 667        ABI_ALLOCATE(veeylm,(ndim,ndim,ndim,ndim))
 668        ABI_ALLOCATE(veeylm2,(tndim,tndim,tndim,tndim))
 669        ABI_ALLOCATE(veerotated,(tndim,tndim,tndim,tndim))
 670        ABI_ALLOCATE(fk,(0:lpawu))
 671        fk(0)=hu(itypat)%upawu
 672        if (lpawu==1) then
 673          fk(1)=hu(itypat)%jpawu*5
 674        else if (lpawu==2) then
 675          fk(1)=hu(itypat)%jpawu*14._dp/(One+hu(itypat)%f4of2_sla)
 676          fk(2)=fk(1)*hu(itypat)%f4of2_sla
 677        else if (lpawu==3) then
 678          fk(1)=hu(itypat)%jpawu*6435._dp/(286._dp+195._dp*hu(itypat)%f4of2_sla+250._dp*hu(itypat)%f6of2_sla)
 679          fk(2)=fk(1)*hu(itypat)%f4of2_sla
 680          fk(3)=fk(1)*hu(itypat)%f6of2_sla
 681        endif
 682 
 683 
 684 
 685 
 686 !      ==================================
 687 !      First define veeslm
 688 !      ==================================
 689 
 690 
 691 !!     veeslm2(s1m1,s2m2,s3m3,s4m4)=vee(m1,m2,m3,m4)*delta_s1s3*delta_s2s4
 692        call vee_ndim2tndim_hu(lpawu,hu(itypat)%vee,veeslm2,1)
 693 
 694 !!     veeslm(m1,m2,m3,m4)=cmplx(vee(m1,m2,m3,m4),zero)
 695        veeslm(:,:,:,:)=cmplx(hu(itypat)%vee(:,:,:,:),zero)
 696 
 697 !!     build udens in the Slm basis and print it
 698        call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,real(veeslm),"slm")
 699 
 700        dim_vee=ndim
 701        basis_vee='Slm'
 702 !     first print veeslm
 703        !call printvee_hu(ndim,real(veeslm),1,basis_vee)
 704        call printvee_hu(tndim,real(veeslm2),1,basis_vee)
 705 !      ==================================
 706 !      Then compute veerotated
 707 !      ==================================
 708 
 709 !      In the basis where levels/densitymatrix/green function is  diagonalized
 710 !      ================================================================================
 711        if (rot_type==1) then
 712 !      ---------------------
 713 
 714          veerotated=czero
 715          do m1=1,tndim
 716            do m2=1,tndim
 717              do m3=1,tndim
 718                do m4=1,tndim
 719                  do mi=1,tndim
 720                    do mj=1,tndim
 721                      do mk=1,tndim
 722                        do ml=1,tndim
 723                           veerotated(m1,m2,m3,m4)= veerotated(m1,m2,m3,m4) + &
 724 &                          conjg(rot_mat(iatom,1)%value(mi,m1))* &
 725 &                          conjg(rot_mat(iatom,1)%value(mj,m2))* &
 726 &                                rot_mat(iatom,1)%value(mk,m3)* &
 727 &                                rot_mat(iatom,1)%value(ml,m4)* &
 728 &                              veeslm2(mi,mj,mk,ml)
 729                        enddo
 730                      enddo
 731                    enddo
 732                  enddo
 733                enddo
 734              enddo
 735            enddo
 736          enddo
 737 
 738          dim_vee=tndim
 739          basis_vee='Diagonal basis from Slm basis'
 740 !      In the Ylm basis
 741 !      ================================================================================
 742        else if (rot_type==2.or.rot_type==3.or.rot_type==4) then
 743 !      ---------------------------
 744 
 745 
 746 !      Change basis from slm to ylm basis
 747          call vee_slm2ylm_hu(lpawu,veeslm,veeylm,1,2)
 748 
 749          basis_vee='Ylm'
 750          dim_vee=ndim
 751 
 752 !        print interaction matrix in the ylm basis
 753          call printvee_hu(ndim,real(veeylm),1,basis_vee,hu(itypat)%upawu)
 754 
 755 !        print interaction matrix in the ylm basis from Slater tables
 756          if(pawprtvol>=3) call udens_slatercondon_hu(fk,lpawu)
 757 
 758          if (rot_type==3.or.rot_type==4) then
 759 !        ---------------------------
 760 !
 761 
 762 !          built large matrix
 763            call vee_ndim2tndim_hu(lpawu,real(veeylm),veeylm2,1)
 764 
 765 
 766            call printvee_hu(tndim,real(veeylm2),1,basis_vee)
 767 
 768          endif
 769 
 770 
 771 !      In the JmJ basis
 772 !      ================================================================================
 773          if (rot_type==3) then
 774 
 775 !          apply change of basis
 776            call vee_ylm2jmj_hu(lpawu,veeylm2,veejmj,1)
 777 
 778 !        print interaction matrix in the JMJ basis from Inglis and Julien tables
 779            if(pawprtvol>=3) call udens_inglis_hu(fk,lpawu)
 780 
 781 !          new dimension
 782 
 783 
 784 
 785            dim_vee=tndim
 786            basis_vee='JmJ'
 787 
 788          else if (rot_type==4) then
 789 
 790            call udens_inglis_hu(fk,lpawu)
 791            veerotated=czero
 792 
 793            do m1=1,tndim
 794              do m2=1,tndim
 795                do m3=1,tndim
 796                  do m4=1,tndim
 797                    do mi=1,tndim
 798                      do mj=1,tndim
 799                        do mk=1,tndim
 800                          do ml=1,tndim
 801                             veerotated(m1,m2,m3,m4)= veerotated(m1,m2,m3,m4) + &
 802 &                          conjg(rot_mat(iatom,1)%value(mi,m1))* &
 803 &                          conjg(rot_mat(iatom,1)%value(mj,m2))* &
 804 &                                rot_mat(iatom,1)%value(mk,m3)* &
 805 &                                rot_mat(iatom,1)%value(ml,m4)* &
 806 &                                veeylm2(mi,mj,mk,ml)
 807                          enddo
 808                        enddo
 809                      enddo
 810                    enddo
 811                  enddo
 812                enddo
 813              enddo
 814            enddo
 815 
 816            dim_vee=tndim
 817            basis_vee='Diagonal basis from Ylm'
 818 
 819          endif
 820        endif
 821        ABI_ALLOCATE(veenew,(dim_vee,dim_vee,dim_vee,dim_vee))
 822        if(rot_type==0) then   ; veenew=veeslm
 823        else if(rot_type==1) then  ; veenew=veerotated
 824        else if(rot_type==2) then  ; veenew=veeylm
 825        else if(rot_type==3) then  ; veenew=veejmj
 826        else if(rot_type==4) then  ; veenew=veerotated
 827        endif
 828 
 829        call printvee_hu(dim_vee,real(veenew),1,basis_vee,hu(itypat)%upawu,hu(itypat)%f2_sla)
 830 !       call printvee_hu(dim_vee,real(veeylm),1,hu(itypat)%upawu)
 831 
 832        uaver=zero
 833        do ms=1,2*ndim
 834          do ms1=1,2*ndim
 835            udens_atoms(iatom)%value(ms,ms1)=veenew(ms,ms1,ms,ms1)-veenew(ms,ms1,ms1,ms)
 836            uaver=uaver+udens_atoms(iatom)%value(ms,ms1)
 837          enddo
 838        enddo
 839        !write(6,*) "u average",uaver/float(4*ndim*ndim)
 840 
 841        call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,real(veenew),"basis_vee",prtonly=1)
 842 
 843 
 844        ABI_DEALLOCATE(veenew)
 845        ABI_DEALLOCATE(veeslm)
 846        ABI_DEALLOCATE(veeslm2)
 847        ABI_DEALLOCATE(veeylm)
 848        ABI_DEALLOCATE(veeylm2)
 849        ABI_DEALLOCATE(veejmj)
 850        ABI_DEALLOCATE(veerotated)
 851        ABI_DEALLOCATE(veetmp_slm)
 852        ABI_DEALLOCATE(veetmp)
 853        ABI_DEALLOCATE(veetmp_ylm)
 854        ABI_DEALLOCATE(fk)
 855      endif
 856    enddo
 857    !MSG_ERROR("Aborting now!")
 858 
 859 !================================================
 860 !  NSPINOR = 1
 861 !================================================
 862  else if (nspinor==1) then
 863 
 864    nat_correl=0
 865    do iatom=1,natom
 866      itypat=cryst_struc%typat(iatom)
 867      lpawu=hu(itypat)%lpawu
 868      if(lpawu.ne.-1) then
 869        nat_correl=nat_correl+1
 870        if(nat_correl>1.and.(hu(itypat)%jpawu>tol4)) then
 871           write(message,'(3a)')  ch10,'  -------> Warning: several atoms: '&
 872 &         ,' not extensively tested '
 873           MSG_WARNING(message)
 874        endif
 875 
 876 !  ! ================================================================
 877 !  !  If rotation for spin 2 and rotation for spin 1 are not equal print
 878 !  !  then print a warning
 879 !  !  useful only for magnetic case
 880 !  ! ================================================================
 881        ndim=2*lpawu+1
 882        tndim=nspinor*ndim
 883        do m1=1,tndim
 884          do m2=1,tndim
 885            if(nsppol==2) then
 886              if(abs(rot_mat(iatom,1)%value(m1,m2)-rot_mat(iatom,2)%value(m1,m2))>tol4.and.pawprtvol>=3) then
 887                write(message,'(2a,i4)')  ch10,' rot_mat differs for value of isppol but value for isppol=2 not used'
 888                call wrtout(std_out,  message,'COLL')
 889                write(message,'(a,4e16.8)')  ch10,rot_mat(iatom,1)%value(m1,m2),rot_mat(iatom,2)%value(m1,m2)
 890                call wrtout(std_out,  message,'COLL', do_flush=.True.)
 891              endif
 892            endif
 893          end do
 894        end do
 895        ABI_ALLOCATE(temp_mat,(ndim,ndim))
 896        ABI_ALLOCATE(temp_mat2,(ndim,ndim))
 897        ABI_ALLOCATE(veetemp,(ndim,ndim,ndim,ndim))
 898        temp_mat(:,:)=czero
 899        temp_mat2(:,:)=czero
 900 
 901 !  ! =================================================
 902 !  !    See if rotation is complex or real
 903 !  ! =================================================
 904        do mi=1,ndim
 905          do m1=1,ndim
 906          if(abs(aimag(rot_mat(iatom,1)%value(mi,m1)))>tol8.and.pawprtvol>=3) then
 907               write(message,'(2a,2i6,2e14.3)')  ch10,"rot_mat is complex for", &
 908 &              mi,m1,rot_mat(iatom,1)%value(mi,m1)
 909               call wrtout(std_out,  message,'COLL')
 910            endif
 911          enddo
 912        enddo
 913 
 914 
 915 
 916 !  !    write vee for information with a classification.
 917        if(pawprtvol>=3) then
 918          call printvee_hu(ndim,hu(itypat)%vee,2,'original')
 919        endif
 920 
 921 !  !    Compute rotated vee.
 922        veetemp=zero
 923        do m1=1,ndim
 924          do m2=1,ndim
 925            do m3=1,ndim
 926              do m4=1,ndim
 927                do mi=1,ndim
 928                  do mj=1,ndim
 929                    do mk=1,ndim
 930                      do ml=1,ndim
 931 !                        if((mi==mk.and.mj==ml).or.(mi==ml.and.mj==mk)) then
 932                         veetemp(m1,m2,m3,m4)= veetemp(m1,m2,m3,m4) + &
 933 &                      real(   &
 934 &                          conjg(rot_mat(iatom,1)%value(mi,m1))* &
 935 &                          conjg(rot_mat(iatom,1)%value(mj,m2))* &
 936 &                                rot_mat(iatom,1)%value(mk,m3)* &
 937 &                                rot_mat(iatom,1)%value(ml,m4)* &
 938 &                            hu(itypat)%vee(mi,mj,mk,ml)&
 939                            )
 940 !                        endif
 941                      enddo
 942                    enddo
 943                  enddo
 944                enddo
 945              enddo
 946            enddo
 947          enddo
 948        enddo
 949        ABI_DEALLOCATE(temp_mat)
 950        ABI_DEALLOCATE(temp_mat2)
 951        xsum=zero
 952        xsum2=zero
 953        xsumnew=zero
 954        xsum2new=zero
 955        do m1=1,ndim
 956          do m2=1,ndim
 957            xsum=xsum+hu(itypat)%vee(m1,m2,m1,m2)
 958            xsum2=xsum2+hu(itypat)%vee(m1,m2,m2,m1)
 959            xsumnew=xsumnew+veetemp(m1,m2,m1,m2)
 960            xsum2new=xsum2new+veetemp(m1,m2,m2,m1)
 961          enddo
 962        enddo
 963        if(abs(xsum-xsumnew)>tol5.or.abs(xsum2-xsum2new)>tol5) then
 964          write(message,'(2a)')  ch10," BUG: New interaction after rotation do not respect sum rules"
 965          call wrtout(std_out,  message,'COLL')
 966          write(message,'(2a,2f14.3)')  ch10,' Comparison of \sum_{m1,m3} vee(m1,m3,m1,m3) before and after rotation is',&
 967 &         xsum,xsumnew
 968          call wrtout(std_out,  message,'COLL')
 969          write(message,'(2a,2f14.3)')  ch10,' Comparison of \sum_{m1,m3} vee(m1,m3,m3,m1) before and after rotation is',&
 970 &         xsum2,xsum2new
 971          call wrtout(std_out,  message,'COLL')
 972        endif
 973        if(pawprtvol>=3) then
 974          write(message,'(2a)')  ch10," VEE ROTATED"
 975          call wrtout(std_out,  message,'COLL')
 976          call printvee_hu(tndim,veetemp,2,'diag1')
 977          write(message,'(a)') ch10
 978          call wrtout(std_out,  message,'COLL')
 979        endif
 980 
 981        write(message,'(2a,i4)')  ch10,'  -------> For Correlated Species', itypat
 982        call wrtout(std_out,  message,'COLL')
 983 
 984        call vee2udensatom_hu(ndim,nspinor,udens_atoms(iatom)%value,veetemp,"diag")
 985 
 986 !       udens_atoms(iatom)%value=zero
 987 !       ij=0
 988 !       do ms=1,2*ndim-1
 989 !         do ms1=ms+1,2*ndim
 990 !           ij=ij+1
 991 !           if(ms<=ndim.and.ms1>ndim) then
 992 !             m1 = ms1 - ndim
 993 !             m  = ms
 994 !             hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)
 995 !             udens_atoms(iatom)%value(ms,ms1)= veetemp(m,m1,m,m1)
 996 !             udens_atoms(iatom)%value(ms1,ms)= udens_atoms(iatom)%value(ms,ms1)
 997 !           else if(ms<=ndim.and.ms1<=ndim) then
 998 !             m1 = ms1
 999 !             m  = ms
1000 !             hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m)
1001 !             udens_atoms(iatom)%value(ms,ms1)= hu(itypat)%uqmc(ij)
1002 !             udens_atoms(iatom)%value(ms1,ms)= udens_atoms(iatom)%value(ms,ms1)
1003 !           else
1004 !             m1 = ms1 - ndim
1005 !             m  = ms  - ndim
1006 !             hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m)
1007 !             udens_atoms(iatom)%value(ms,ms1)= hu(itypat)%uqmc(ij)
1008 !             udens_atoms(iatom)%value(ms1,ms)= udens_atoms(iatom)%value(ms,ms1)
1009 !           endif
1010 !         enddo
1011 !       enddo
1012 !       write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation "
1013 !       call wrtout(std_out,  message,'COLL')
1014 !       write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim)
1015 !       call wrtout(std_out,  message,'COLL')
1016 !       do ms=1,2*ndim
1017 !          write(message,'(i3,14f7.3)') &
1018 !  &        ms, (udens_atoms(iatom)%value(ms,ms1),ms1=1,2*ndim)
1019 !          call wrtout(std_out,  message,'COLL')
1020 !       enddo
1021 !       write(message,'(5x,a)') "--------------------------------------------------------"
1022 !       call wrtout(std_out,  message,'COLL')
1023        ABI_DEALLOCATE(veetemp)
1024      endif ! lpawu/=1
1025 !   call print_hu(hu,cryst_struc%ntypat,1)
1026 
1027    enddo ! iatom
1028 !   call print_hu(hu,cryst_struc%ntypat,1)
1029 !   call vee2udens_hu(hu,cryst_struc%ntypat,2)
1030  endif ! nspinor==1
1031 
1032 
1033 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.

PARENTS

CHILDREN

      wrtout

SOURCE

473 subroutine vee2udens_hu(hu,ntypat,prtopt)
474 
475  use defs_basis
476  use m_crystal, only : crystal_t
477 
478 !This section has been created automatically by the script Abilint (TD).
479 !Do not modify the following lines by hand.
480 #undef ABI_FUNC
481 #define ABI_FUNC 'vee2udens_hu'
482 !End of the abilint section
483 
484  implicit none
485 
486 !Arguments ------------------------------------
487 !type
488  integer, intent(in):: ntypat
489  type(hu_type),intent(inout) :: hu(ntypat)
490  integer :: prtopt
491 
492 !Local variables-------------------------------
493  integer :: ij,itypat
494  integer :: lpawu,m1,ms,ms1,m,ndim
495  character(len=500) :: message
496 ! *********************************************************************
497  do itypat=1,ntypat
498    lpawu=hu(itypat)%lpawu
499    if(lpawu.ne.-1) then
500      ndim=2*lpawu+1
501      write(message,'(2a,i4)')  ch10,'  -------> For Correlated Species', itypat
502      call wrtout(std_out,  message,'COLL')
503 
504      hu(itypat)%udens=zero
505      ij=0
506      do ms=1,2*ndim-1
507 !         xij(ms,ms)=0
508        do ms1=ms+1,2*ndim
509          ij=ij+1
510 !         xij(ms,ms1)=ij
511 !         xij(ms1,ms)=ij
512          if(ms<=ndim.and.ms1>ndim) then
513            m1 = ms1 - ndim
514            m  = ms
515            hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)
516            hu(itypat)%udens(ms,ms1)= hu(itypat)%vee(m,m1,m,m1)
517            hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1)
518          else if(ms<=ndim.and.ms1<=ndim) then
519            m1 = ms1
520            m  = ms
521            hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m)
522            hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij)
523            hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1)
524          else
525            m1 = ms1 - ndim
526            m  = ms  - ndim
527            hu(itypat)%uqmc(ij)=hu(itypat)%vee(m,m1,m,m1)-hu(itypat)%vee(m,m1,m1,m)
528            hu(itypat)%udens(ms,ms1)= hu(itypat)%uqmc(ij)
529            hu(itypat)%udens(ms1,ms)= hu(itypat)%udens(ms,ms1)
530          endif
531        enddo
532      enddo
533 !     xij(2*ndim,2*ndim)=0
534 !     write(message,'(a,5x,a)') ch10,"-------- Interactions in the density matrix representation "
535 !     call wrtout(std_out,  message,'COLL')
536 !     write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim)
537 !     call wrtout(std_out,  message,'COLL')
538    endif
539  enddo ! itypat
540  call print_hu(hu,ntypat,prtopt)
541 
542 
543 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

1416 subroutine vee2udensatom_hu(ndim,nspinor,udens_atoms,veetemp,basis,prtonly)
1417 
1418  use defs_basis
1419  use m_crystal, only : crystal_t
1420 
1421 !This section has been created automatically by the script Abilint (TD).
1422 !Do not modify the following lines by hand.
1423 #undef ABI_FUNC
1424 #define ABI_FUNC 'vee2udensatom_hu'
1425 !End of the abilint section
1426 
1427  implicit none
1428 
1429 !Arguments ------------------------------------
1430 !type
1431  integer, intent(in):: ndim,nspinor
1432  real(dp),intent(out) :: udens_atoms(2*ndim,2*ndim)
1433  !real(dp), intent(in) :: veetemp(nspinor*ndim,nspinor*ndim,nspinor*ndim,nspinor*ndim)
1434  real(dp), intent(in) :: veetemp(ndim,ndim,ndim,ndim)
1435  character(len=*), intent(in) :: basis
1436  integer, intent(in), optional:: prtonly
1437 
1438 !Local variables-------------------------------
1439  integer :: ij,ms,ms1,m,m1,prt
1440  character(len=1000) :: message
1441 ! *********************************************************************
1442  prt=1
1443  if(present(prtonly)) prt=2
1444  if (nspinor==1.or.nspinor==2.and.prt<=1) then
1445    udens_atoms=zero
1446    ij=0
1447    do ms=1,2*ndim-1
1448      do ms1=ms+1,2*ndim
1449        ij=ij+1
1450        if(ms<=ndim.and.ms1>ndim) then
1451          m1 = ms1 - ndim
1452          m  = ms
1453 !         hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)
1454          udens_atoms(ms,ms1)= veetemp(m,m1,m,m1)
1455          udens_atoms(ms1,ms)= udens_atoms(ms,ms1)
1456         ! write(6,*)"A", ms,ms1,udens_atoms(ms,ms1)
1457        else if(ms<=ndim.and.ms1<=ndim) then
1458          m1 = ms1
1459          m  = ms
1460 !         hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m)
1461          udens_atoms(ms,ms1)= veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m)
1462          udens_atoms(ms1,ms)= udens_atoms(ms,ms1)
1463         ! write(6,*)"B", ms,ms1,udens_atoms(ms,ms1)
1464        else
1465          m1 = ms1 - ndim
1466          m  = ms  - ndim
1467 !         hu(itypat)%uqmc(ij)=veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m)
1468          udens_atoms(ms,ms1)= veetemp(m,m1,m,m1)-veetemp(m,m1,m1,m)
1469          udens_atoms(ms1,ms)= udens_atoms(ms,ms1)
1470         ! write(6,*)"C", ms,ms1,udens_atoms(ms,ms1)
1471        endif
1472      enddo
1473    enddo
1474 
1475 ! else if(nspinor==2) then
1476 !
1477 !   do ms=1,2*ndim
1478 !     do ms1=1,2*ndim
1479 !       udens_atoms(ms,ms1)=veetemp(ms,ms1,ms,ms1)-veetemp(ms,ms1,ms1,ms)
1480 !     enddo
1481 !   enddo
1482 
1483  endif
1484 
1485 
1486  write(message,'(4a)') ch10,"-------- Interactions in the ",trim(basis)," basis "
1487  call wrtout(std_out,  message,'COLL')
1488  write(message,'(1x,14(2x,i5))') (m,m=1,2*ndim)
1489  call wrtout(std_out,  message,'COLL')
1490  do ms=1,2*ndim
1491     write(message,'(i3,14f7.3)') &
1492 &    ms, (udens_atoms(ms,ms1),ms1=1,2*ndim)
1493     call wrtout(std_out,  message,'COLL')
1494  enddo
1495  write(message,'(5x,a)') "--------------------------------------------------------"
1496  call wrtout(std_out,  message,'COLL')
1497 
1498 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-2018 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

1826 subroutine vee_ndim2tndim_hu(lcor,mat_inp_c,mat_out_c,option)
1827 
1828  use defs_basis
1829  use m_errors
1830  use m_abicore
1831 
1832 !This section has been created automatically by the script Abilint (TD).
1833 !Do not modify the following lines by hand.
1834 #undef ABI_FUNC
1835 #define ABI_FUNC 'vee_ndim2tndim_hu'
1836 !End of the abilint section
1837 
1838  implicit none
1839 
1840 !Arguments ---------------------------------------------
1841 !scalars
1842  integer,intent(in) :: lcor,option
1843 !arrays
1844  real(dp), intent(in) :: mat_inp_c(2*lcor+1,2*lcor+1,2*lcor+1,2*lcor+1)
1845  complex(dpc), intent(out) :: mat_out_c(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1))
1846 
1847 !Local variables ---------------------------------------
1848 !scalars
1849  integer :: m1,m2,m3,m4,is1,is2,is3,is4,ndim,s1,s2,s3,s4
1850 
1851 ! *********************************************************************
1852   ndim=2*lcor+1
1853   mat_out_c=czero
1854   if(option==1) then
1855     do m1=1,ndim
1856       do m2=1,ndim
1857         do m3=1,ndim
1858           do m4=1,ndim
1859             do is1=1,2
1860               do is2=1,2
1861 
1862                 is3=is1 ; is4=is2
1863 
1864                 s1=(is1-1)*ndim ; s2=(is2-1)*ndim ; s3=(is3-1)*ndim ; s4=(is4-1)*ndim
1865 
1866                 mat_out_c(m1+s1,m2+s2,m3+s3,m4+s4)=  cmplx(mat_inp_c(m1,m2,m3,m4),zero)
1867 
1868               enddo
1869             enddo
1870           enddo
1871         enddo
1872       enddo
1873     enddo
1874   endif
1875 
1876 
1877 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-2018 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

PARENTS

      qmc_prep_ctqmc

CHILDREN

      wrtout

SOURCE

1736 subroutine vee_ndim2tndim_hu_r(lcor,mat_inp_c,mat_out_c,option)
1737 
1738  use defs_basis
1739  use m_errors
1740  use m_abicore
1741 
1742 !This section has been created automatically by the script Abilint (TD).
1743 !Do not modify the following lines by hand.
1744 #undef ABI_FUNC
1745 #define ABI_FUNC 'vee_ndim2tndim_hu_r'
1746 !End of the abilint section
1747 
1748  implicit none
1749 
1750 !Arguments ---------------------------------------------
1751 !scalars
1752  integer,intent(in) :: lcor,option
1753 !arrays
1754  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(:,:,:,:)
1755  real(dp), intent(out)       :: mat_out_c(:,:,:,:) !(2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1),2*(2*lcor+1))
1756 
1757 !Local variables ---------------------------------------
1758 !scalars
1759  integer :: m1,m2,m3,m4,is1,is2,is3,is4,ndim,s1,s2,s3,s4
1760 
1761 ! *********************************************************************
1762   ndim=2*lcor+1
1763   mat_out_c=czero
1764 
1765   if(option==1) then
1766     do m1=1,ndim
1767       do m2=1,ndim
1768         do m3=1,ndim
1769           do m4=1,ndim
1770             do is1=1,2
1771               do is2=1,2
1772 
1773                 is3=is1 ; is4=is2
1774 
1775                 s1=(is1-1)*ndim ; s2=(is2-1)*ndim ; s3=(is3-1)*ndim ; s4=(is4-1)*ndim
1776 
1777                 mat_out_c(m1+s1,m2+s2,m3+s3,m4+s4)=  mat_inp_c(m1,m2,m3,m4)
1778 
1779               enddo
1780             enddo
1781           enddo
1782         enddo
1783       enddo
1784     enddo
1785   endif
1786 
1787 
1788 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-2018 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

PARENTS

      m_hu

CHILDREN

      wrtout

SOURCE

1578 subroutine vee_slm2ylm_hu(lcor,mat_inp_c,mat_out_c,option,prtvol)
1579 
1580  use defs_basis
1581  use m_errors
1582  use m_abicore
1583 
1584 !This section has been created automatically by the script Abilint (TD).
1585 !Do not modify the following lines by hand.
1586 #undef ABI_FUNC
1587 #define ABI_FUNC 'vee_slm2ylm_hu'
1588 !End of the abilint section
1589 
1590  implicit none
1591 
1592 !Arguments ---------------------------------------------
1593 !scalars
1594  integer,intent(in) :: lcor,option,prtvol
1595 !arrays
1596  complex(dpc), intent(in) :: mat_inp_c(2*lcor+1,2*lcor+1,2*lcor+1,2*lcor+1)
1597  complex(dpc), intent(out) :: mat_out_c(2*lcor+1,2*lcor+1,2*lcor+1,2*lcor+1)
1598 
1599 !Local variables ---------------------------------------
1600 !scalars
1601  integer :: gm,hm,jm,gg,hh,ii,jj,ll,mm,im
1602  real(dp),parameter :: invsqrt2=one/sqrt2
1603  real(dp) :: onem
1604  complex(dpc) :: tmp2
1605  character(len=500) :: message
1606 !arrays
1607  complex(dpc),allocatable :: slm2ylm(:,:)
1608 
1609 ! *********************************************************************
1610 
1611 
1612  if (option/=1.and.option/=2) then
1613    message=' option=/1 or 2 !'
1614    MSG_BUG(message)
1615  end if
1616 
1617  if(abs(prtvol)>2) then
1618    write(message,'(3a)') ch10, "   vee_slm2ylm_hu"
1619    call wrtout(std_out,message,'COLL')
1620  end if
1621 
1622  if(abs(prtvol)>2) then
1623    if(option==1) then
1624      write(message,'(3a)') ch10,"matrix in Slm basis is changed into Ylm basis"
1625      call wrtout(std_out,message,'COLL')
1626    else if(option==2) then
1627      write(message,'(3a)') ch10,"matrix in Ylm basis is changed into Slm basis"
1628      call wrtout(std_out,message,'COLL')
1629    end if
1630  end if
1631 
1632  ll=lcor
1633  ABI_ALLOCATE(slm2ylm,(2*ll+1,2*ll+1))
1634  slm2ylm=czero
1635  mat_out_c=czero
1636 
1637 !  ===== Definitions of slm2ylm
1638  do im=1,2*ll+1
1639    mm=im-ll-1;jm=-mm+ll+1   ! mmj=-mm
1640 !! im is in {1,....2*ll+1}
1641 !! mm is in {-ll,....+ll}
1642 !! jm is in {2*ll+1,....,1}
1643    onem=dble((-1)**mm)
1644    if (mm> 0) then ! im in {ll+1,2ll+1} and jm in {ll+1,1}
1645      slm2ylm(im,im)= cmplx(onem*invsqrt2,zero,kind=dp)
1646      slm2ylm(jm,im)= cmplx(invsqrt2,     zero,kind=dp)
1647    end if
1648    if (mm==0) then
1649      slm2ylm(im,im)=cone
1650    end if
1651    if (mm< 0) then
1652      slm2ylm(im,im)= cmplx(zero,     invsqrt2,kind=dp)
1653      slm2ylm(jm,im)=-cmplx(zero,onem*invsqrt2,kind=dp)
1654    end if
1655  end do
1656 ! do im=1,2*ll+1
1657 !   write(message,'(7(2f14.5))') (slm2ylm(im,jm),jm=1,2*ll+1)
1658 !   call wrtout(std_out,message,'COLL')
1659 ! end do
1660 
1661 !  ===== Definitions of slm2ylm
1662 !!!!  pawtab(itypat)%vee(m11,m31,m21,m41)= <m11 m31| vee| m21 m41 >
1663 !!!!  pawtab(itypat)%vee(m11,m21,m31,m41)= <m11 m21| vee| m31 m41 >
1664 
1665  do jm=1,2*ll+1
1666    do im=1,2*ll+1
1667      do hm=1,2*ll+1
1668        do gm=1,2*ll+1
1669          tmp2=czero
1670          do gg=1,2*ll+1
1671            do hh=1,2*ll+1
1672              do ii=1,2*ll+1
1673                do jj=1,2*ll+1
1674                  if(option==1) then
1675                    tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*(slm2ylm(im,ii))*CONJG(slm2ylm(jm,jj))&
1676 &                                                   *(slm2ylm(gm,gg))*CONJG(slm2ylm(hm,hh))
1677 !                   if(gm==1.and.hm==1.and.im==1.and.jm==1) then
1678 !                      write(6,'(4i4,2f10.5,2f10.5)') gg,hh,ii,jj,tmp2,mat_inp_c(gg,hh,ii,jj)
1679 !                      write(6,*) "i1"
1680 !                   endif
1681                  else if(option==2) then
1682                    tmp2=tmp2+mat_inp_c(gg,ii,hh,jj)*CONJG(slm2ylm(ii,im))*(slm2ylm(jj,jm))&
1683 &                                                   *CONJG(slm2ylm(gg,gm))*(slm2ylm(hh,hm))
1684                  end if
1685                end do
1686              end do
1687            end do
1688          end do
1689 !         mat_out_c(gm,hm,im,jm)=tmp2
1690          mat_out_c(gm,im,hm,jm)=tmp2
1691        end do
1692      end do
1693    end do
1694  end do
1695 
1696  ABI_DEALLOCATE(slm2ylm)
1697 
1698 end subroutine vee_slm2ylm_hu