TABLE OF CONTENTS
- ABINIT/m_plowannier
- m_plowannier/allocate_orbital
- m_plowannier/compute_oper_ks2wan
- m_plowannier/copy_orbital
- m_plowannier/destroy_operwan
- m_plowannier/destroy_orbital
- m_plowannier/destroy_plowannier
- m_plowannier/initialize_operwan
- m_plowannier/latom_wan_type
- m_plowannier/lorbital2_type
- m_plowannier/lorbital_type
- m_plowannier/normalization_plowannier
- m_plowannier/operwan_realspace_type
- m_plowannier/operwan_type
- m_plowannier/orbital_type
- m_plowannier/plowannier_type
- m_plowannier/position_wan_type
- m_plowannier/print_operwan
- m_plowannier/projector_wan_type
ABINIT/m_plowannier [ Modules ]
NAME
m_plowannier
FUNCTION
COPYRIGHT
Copyright (C) 2006-2018 ABINIT group (BAmadon,AGerossier) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
PARENTS
NOTES
SOURCE
20 #if defined HAVE_CONFIG_H 21 #include "config.h" 22 #endif 23 24 #include "abi_common.h" 25 26 MODULE m_plowannier 27 28 use defs_basis 29 use m_errors 30 use m_profiling_abi 31 32 use m_io_tools, only : open_file 33 34 implicit none 35 36 private 37 38 public :: init_plowannier 39 public :: copy_orbital 40 public :: compute_coeff_plowannier 41 public :: destroy_plowannier 42 public :: initialize_operwan 43 public :: destroy_operwan 44 public :: compute_oper_ks2wan 45 public :: normalization_plowannier 46 public :: print_operwan
m_plowannier/allocate_orbital [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
allocate_orbital
FUNCTION
allocate an array of orbital_type
INPUTS
lorbital1
OUTPUT
lorbital2
PARENTS
m_plowannier
CHILDREN
SOURCE
552 subroutine allocate_orbital(orbital1,orbital2,n1,n2,n3) 553 554 555 !This section has been created automatically by the script Abilint (TD). 556 !Do not modify the following lines by hand. 557 #undef ABI_FUNC 558 #define ABI_FUNC 'allocate_orbital' 559 !End of the abilint section 560 561 implicit none 562 563 !Arguments---------------- 564 integer,intent(in) :: n1,n2,n3 565 type(orbital_type), intent(in) :: orbital1(n1,n2,n3) 566 type(orbital_type),intent(inout) :: orbital2(n1,n2,n3) 567 568 !Local variable----------- 569 integer :: n4,n5,n6,n7 570 integer :: i,j,k,l 571 572 do i = 1,n1 573 do j = 1,n2 574 do k = 1,n3 575 n4 = size(orbital1(i,j,k)%atom,1) 576 ABI_DATATYPE_ALLOCATE(orbital2(i,j,k)%atom,(n4)) 577 do l = 1,n4 578 n5 = size(orbital1(i,j,k)%atom(l)%matl,1) 579 n6 = size(orbital1(i,j,k)%atom(l)%matl,2) 580 n7 = size(orbital1(i,j,k)%atom(l)%matl,3) 581 ABI_ALLOCATE(orbital2(i,j,k)%atom(l)%matl,(n5,n6,n7)) 582 end do 583 end do 584 end do 585 end do 586 587 end subroutine allocate_orbital
m_plowannier/compute_oper_ks2wan [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
compute_oper_ks2wan
FUNCTION
transform ks operator into wan one
INPUTS
wan,operks,option
OUTPUT
if option = ikpt, gives the wan operator in reciprocal space (for each k)
PARENTS
m_plowannier
CHILDREN
SOURCE
2411 subroutine compute_oper_ks2wan(wan,operks,operwan,option) 2412 2413 #ifdef FC_INTEL 2414 !DEC$ NOOPTIMIZE 2415 #endif 2416 2417 !This section has been created automatically by the script Abilint (TD). 2418 !Do not modify the following lines by hand. 2419 #undef ABI_FUNC 2420 #define ABI_FUNC 'compute_oper_ks2wan' 2421 !End of the abilint section 2422 2423 implicit none 2424 2425 !Arguments-------------------------- 2426 type(plowannier_type), intent(in) :: wan 2427 type(operwan_type), intent(inout) :: operwan(:,:,:) 2428 complex(dpc), intent(in) :: operks(:,:,:,:) 2429 integer, intent(in) :: option 2430 2431 !Local variables-------------------- 2432 integer :: iatom1, iatom2, il1, il2, isppol, ispinor1, ispinor2, iband1, iband2, im1, im2 2433 2434 ! ---------------------------------- 2435 !Transformation KS2WAN 2436 ! ---------------------------------- 2437 2438 2439 !!operation on reciprocal space 2440 do iatom1 = 1,wan%natom_wan 2441 do iatom2 = 1,wan%natom_wan 2442 do isppol = 1,wan%nsppol 2443 do il1 = 1,wan%nbl_atom_wan(iatom1) 2444 do il2 = 1,wan%nbl_atom_wan(iatom2) 2445 do im1 = 1,2*wan%latom_wan(iatom1)%lcalc(il1)+1 2446 do im2 = 1,2*wan%latom_wan(iatom2)%lcalc(il2)+1 2447 do ispinor1 = 1,wan%nspinor 2448 do ispinor2 = 1,wan%nspinor 2449 !!sum over the bands 2450 do iband1 = 1,wan%bandf_wan-wan%bandi_wan+1 2451 do iband2 = 1,wan%bandf_wan-wan%bandi_wan+1 2452 operwan(option,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2)=& 2453 operwan(option,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2)& 2454 + conjg(wan%psichi(option,iband2,iatom2)%atom(il2)%matl(im2,isppol,ispinor2))& 2455 *operks(option,iband1,iband2,isppol)*wan%psichi(option,iband1,iatom1)%atom(il1)%matl(im1,isppol,ispinor1) 2456 end do 2457 end do 2458 end do 2459 end do 2460 end do 2461 end do 2462 end do 2463 end do 2464 end do 2465 end do 2466 end do 2467 2468 end subroutine compute_oper_ks2wan
m_plowannier/copy_orbital [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
copy_orbital
FUNCTION
Copy an array of orbital_type
INPUTS
lorbital1
OUTPUT
lorbital2
PARENTS
m_plowannier
CHILDREN
SOURCE
487 subroutine copy_orbital(orbital1,orbital2,n1,n2,n3) 488 489 490 !This section has been created automatically by the script Abilint (TD). 491 !Do not modify the following lines by hand. 492 #undef ABI_FUNC 493 #define ABI_FUNC 'copy_orbital' 494 !End of the abilint section 495 496 implicit none 497 498 !Arguments---------------- 499 integer,intent(in) :: n1,n2,n3 500 type(orbital_type), intent(in) :: orbital1(n1,n2,n3) 501 type(orbital_type),intent(inout) :: orbital2(n1,n2,n3) 502 503 !Local variable----------- 504 integer :: n4,n5,n6,n7 505 integer :: i,j,k,l,m,p,q 506 507 do i = 1,n1 508 do j = 1,n2 509 do k = 1,n3 510 n4 = size(orbital1(i,j,k)%atom,1) 511 do l = 1,n4 512 n5 = size(orbital1(i,j,k)%atom(l)%matl,1) 513 n6 = size(orbital1(i,j,k)%atom(l)%matl,2) 514 n7 = size(orbital1(i,j,k)%atom(l)%matl,3) 515 do m = 1,n5 516 do p = 1,n6 517 do q = 1,n7 518 orbital2(i,j,k)%atom(l)%matl(m,p,q) = orbital1(i,j,k)%atom(l)%matl(m,p,q) 519 end do 520 end do 521 end do 522 end do 523 end do 524 end do 525 end do 526 527 end subroutine copy_orbital
m_plowannier/destroy_operwan [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
destroy_operwan
FUNCTION
destroy operwan
INPUTS
wan
OUTPUT
operwan
PARENTS
m_plowannier
CHILDREN
SOURCE
2356 subroutine destroy_operwan(wan,operwan) 2357 2358 2359 !This section has been created automatically by the script Abilint (TD). 2360 !Do not modify the following lines by hand. 2361 #undef ABI_FUNC 2362 #define ABI_FUNC 'destroy_operwan' 2363 !End of the abilint section 2364 2365 implicit none 2366 2367 !Arguments---------------------------------- 2368 type(plowannier_type), intent(in) :: wan 2369 type(operwan_type), intent(inout) :: operwan(wan%nkpt,wan%natom_wan,wan%natom_wan) 2370 2371 !Local variables---------------------------- 2372 integer :: ikpt,iatom1,iatom2,il1,il2 2373 2374 2375 do ikpt = 1,wan%nkpt 2376 do iatom1 = 1,wan%natom_wan 2377 do iatom2 = 1,wan%natom_wan 2378 do il1 = 1,wan%nbl_atom_wan(iatom1) 2379 do il2 = 1,wan%nbl_atom_wan(iatom2) 2380 ABI_DEALLOCATE(operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl) 2381 end do 2382 end do 2383 ABI_DATATYPE_DEALLOCATE(operwan(ikpt,iatom1,iatom2)%atom) 2384 end do 2385 end do 2386 end do 2387 end subroutine destroy_operwan
m_plowannier/destroy_orbital [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
destroy_orbital
FUNCTION
destroy an array of orbital_type
INPUTS
lorbital1
OUTPUT
lorbital2
PARENTS
m_plowannier
CHILDREN
SOURCE
613 subroutine destroy_orbital(orbital2,n1,n2,n3) 614 615 616 !This section has been created automatically by the script Abilint (TD). 617 !Do not modify the following lines by hand. 618 #undef ABI_FUNC 619 #define ABI_FUNC 'destroy_orbital' 620 !End of the abilint section 621 622 implicit none 623 624 !Arguments---------------- 625 integer,intent(in) :: n1,n2,n3 626 type(orbital_type),intent(inout) :: orbital2(n1,n2,n3) 627 628 !Local variable----------- 629 integer :: n4 630 integer :: i,j,k,l 631 632 do i = 1,n1 633 do j = 1,n2 634 do k = 1,n3 635 n4 = size(orbital2(i,j,k)%atom,1) 636 do l = 1,n4 637 ABI_DEALLOCATE(orbital2(i,j,k)%atom(l)%matl) 638 end do 639 ABI_DATATYPE_DEALLOCATE(orbital2(i,j,k)%atom) 640 end do 641 end do 642 end do 643 644 end subroutine destroy_orbital
m_plowannier/destroy_plowannier [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
destroy_plowannier
FUNCTION
deallocate variables
INPUTS
wan
OUTPUT
PARENTS
outscfcv
CHILDREN
SOURCE
2199 subroutine destroy_plowannier(wan) 2200 2201 2202 !This section has been created automatically by the script Abilint (TD). 2203 !Do not modify the following lines by hand. 2204 #undef ABI_FUNC 2205 #define ABI_FUNC 'destroy_plowannier' 2206 !End of the abilint section 2207 2208 implicit none 2209 2210 !Arguments------------------------------------- 2211 type(plowannier_type), intent(inout) :: wan 2212 !Local variables------------------------------- 2213 integer :: iatom,ikpt,iband,il 2214 2215 do iatom=1,wan%natom_wan 2216 ABI_DEALLOCATE(wan%latom_wan(iatom)%lcalc) 2217 ABI_DEALLOCATE(wan%projector_wan(iatom)%lproj) 2218 ABI_DEALLOCATE(wan%nposition(iatom)%pos) 2219 enddo 2220 do iatom = 1,wan%natom_wan 2221 do il = 1,wan%nbl_atom_wan(iatom) 2222 ABI_DEALLOCATE(wan%psichi(1,1,iatom)%atom(il)%ph0phiint) 2223 end do 2224 end do 2225 do ikpt = 1,wan%nkpt 2226 do iband = wan%bandi_wan,wan%bandf_wan 2227 do iatom = 1,wan%natom_wan 2228 do il = 1,wan%nbl_atom_wan(iatom) 2229 ABI_DEALLOCATE(wan%psichi(ikpt,iband-wan%bandi_wan+1,iatom)%atom(il)%matl) 2230 end do 2231 ABI_DATATYPE_DEALLOCATE(wan%psichi(ikpt,iband-wan%bandi_wan+1,iatom)%atom) 2232 end do 2233 end do 2234 end do 2235 2236 if (allocated(wan%kpt)) then 2237 ABI_DEALLOCATE(wan%kpt) 2238 end if 2239 if (allocated(wan%iatom_wan)) then 2240 ABI_DEALLOCATE(wan%iatom_wan) 2241 end if 2242 if (allocated(wan%nbl_atom_wan)) then 2243 ABI_DEALLOCATE(wan%nbl_atom_wan) 2244 end if 2245 if (allocated(wan%latom_wan)) then 2246 ABI_DATATYPE_DEALLOCATE(wan%latom_wan) 2247 end if 2248 if (allocated(wan%nbproj_atom_wan)) then 2249 ABI_DEALLOCATE(wan%nbproj_atom_wan) 2250 end if 2251 if (allocated(wan%projector_wan)) then 2252 ABI_DATATYPE_DEALLOCATE(wan%projector_wan) 2253 end if 2254 if (allocated(wan%position)) then 2255 ABI_DEALLOCATE(wan%position) 2256 end if 2257 if (allocated(wan%wtk)) then 2258 ABI_DEALLOCATE(wan%wtk) 2259 end if 2260 if (allocated(wan%acell)) then 2261 ABI_DEALLOCATE(wan%acell) 2262 end if 2263 2264 2265 if (allocated(wan%nposition)) then 2266 ABI_DATATYPE_DEALLOCATE(wan%nposition) 2267 end if 2268 if (allocated(wan%psichi)) then 2269 ABI_DATATYPE_DEALLOCATE(wan%psichi) 2270 end if 2271 2272 2273 end subroutine destroy_plowannier
m_plowannier/initialize_operwan [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
initialize_operwan
FUNCTION
initialize operwan
INPUTS
wan
OUTPUT
operwan
PARENTS
m_plowannier
CHILDREN
SOURCE
2298 subroutine initialize_operwan(wan,operwan) 2299 2300 2301 !This section has been created automatically by the script Abilint (TD). 2302 !Do not modify the following lines by hand. 2303 #undef ABI_FUNC 2304 #define ABI_FUNC 'initialize_operwan' 2305 !End of the abilint section 2306 2307 implicit none 2308 2309 !Arguments---------------------------------- 2310 type(plowannier_type), intent(in) :: wan 2311 type(operwan_type), intent(inout) :: operwan(wan%nkpt,wan%natom_wan,wan%natom_wan) 2312 2313 !Local variables---------------------------- 2314 integer :: ikpt,iatom1,iatom2,il1,il2,n1,n2 2315 2316 do ikpt = 1,wan%nkpt 2317 do iatom1 = 1,wan%natom_wan 2318 do iatom2 = 1,wan%natom_wan 2319 ABI_DATATYPE_ALLOCATE(operwan(ikpt,iatom1,iatom2)%atom,(wan%nbl_atom_wan(iatom1),wan%nbl_atom_wan(iatom2))) 2320 do il1 = 1,wan%nbl_atom_wan(iatom1) 2321 do il2 = 1,wan%nbl_atom_wan(iatom2) 2322 n1=2*wan%latom_wan(iatom1)%lcalc(il1)+1 2323 n2=2*wan%latom_wan(iatom2)%lcalc(il2)+1 2324 ABI_ALLOCATE(operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl,(n1,n2,wan%nsppol,wan%nspinor,wan%nspinor)) 2325 operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl = zero 2326 end do 2327 end do 2328 end do 2329 end do 2330 end do 2331 2332 end subroutine initialize_operwan
m_plowannier/latom_wan_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
latom_wan_type
FUNCTION
SOURCE
59 type, public :: latom_wan_type 60 61 integer, allocatable :: lcalc(:) 62 ! array of the l we want to compute the psichi with 63 64 end type latom_wan_type
m_plowannier/lorbital2_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
lorbital2_type
FUNCTION
SOURCE
152 type, public :: lorbital2_type 153 154 complex(dpc), allocatable :: matl(:,:,:,:,:) 155 ! size (2l1+1,2l2+1,nspppol,nspinor,nspinor) 156 157 real(dp), allocatable :: ph0phiint(:) 158 ! size (nproj), stocks the value of ph0phiint we may want 159 160 161 end type lorbital2_type
m_plowannier/lorbital_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
lorbital_type
FUNCTION
SOURCE
112 type, public :: lorbital_type 113 114 complex(dpc), allocatable :: matl(:,:,:) 115 !details for different m 116 117 real(dp), allocatable :: ph0phiint(:) 118 ! stocks the values for each projector of the l considered 119 ! size(total number of projectors for this l) 120 121 end type lorbital_type
m_plowannier/normalization_plowannier [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
normalization_plowannier
FUNCTION
Use compute_oper_ks2wan to calculate overlap and do the normalization for the wan%psichi coefficients
INPUTS
wan, option (=0 if we normalize over all k)
OUTPUT
wan itself is modified
PARENTS
m_plowannier
CHILDREN
SOURCE
2494 subroutine normalization_plowannier(wan) 2495 2496 2497 use m_matrix, only : invsqrt_matrix 2498 2499 !This section has been created automatically by the script Abilint (TD). 2500 !Do not modify the following lines by hand. 2501 #undef ABI_FUNC 2502 #define ABI_FUNC 'normalization_plowannier' 2503 !End of the abilint section 2504 2505 implicit none 2506 !Arguments------------------ 2507 2508 !This section has been created automatically by the script Abilint (TD). 2509 !Do not modify the following lines by hand. 2510 #undef ABI_FUNC 2511 #define ABI_FUNC 'normalization_plowannier' 2512 !End of the abilint section 2513 2514 type(plowannier_type), intent(inout) :: wan 2515 2516 !Local---------------------- 2517 complex(dpc), allocatable :: operks(:,:,:,:) 2518 type(operwan_type), allocatable :: operwan(:,:,:) 2519 complex(dpc), allocatable :: operwansquare(:,:,:,:) 2520 complex(dpc), allocatable :: tmp_operwansquare(:,:) 2521 integer :: ikpt, iband, iband1, iband2, isppol, ispinor1, ispinor2, iatom1, & 2522 & iatom2, il1, il2, im1, im2, index_c, index_l, n1,n2,n3 2523 type(orbital_type), allocatable :: psichinormalized(:,:,:) 2524 !character(len = 50) :: mat_writing2 2525 !character(len = 5000) :: mat_writing 2526 character(len = 500) :: message 2527 2528 2529 !First, creation of the ks identity operator 2530 ABI_ALLOCATE(operks,(wan%nkpt,wan%bandf_wan-wan%bandi_wan+1,wan%bandf_wan-wan%bandi_wan+1,wan%nsppol)) 2531 operks = czero 2532 do iband1 = 1,wan%bandf_wan-wan%bandi_wan+1 2533 do iband2 = 1,wan%bandf_wan-wan%bandi_wan+1 2534 if (iband1.eq.iband2) then 2535 do ikpt = 1,wan%nkpt 2536 do isppol= 1,wan%nsppol 2537 operks(ikpt,iband1,iband2,isppol) = cone 2538 end do 2539 end do 2540 end if 2541 end do 2542 end do 2543 2544 2545 !Allocation of operwan 2546 ABI_DATATYPE_ALLOCATE(operwan,(wan%nkpt,wan%natom_wan,wan%natom_wan)) 2547 call initialize_operwan(wan,operwan) 2548 2549 2550 2551 !Computation of the overlap 2552 do ikpt = 1,wan%nkpt 2553 call compute_oper_ks2wan(wan,operks,operwan,ikpt) 2554 end do 2555 2556 2557 2558 !transform the operwan into an inversible matrix 2559 !!operwansquare is the overlap square matrix (wan%size_wan * wan%size_wan) 2560 ABI_ALLOCATE(operwansquare,(wan%nkpt,wan%nsppol,wan%nspinor*wan%size_wan,wan%nspinor*wan%size_wan)) 2561 2562 operwansquare = czero 2563 2564 n1=size(wan%psichi,1) 2565 n2=size(wan%psichi,2) 2566 n3=size(wan%psichi,3) 2567 ABI_DATATYPE_ALLOCATE(psichinormalized,(n1,n2,n3)) 2568 call allocate_orbital(wan%psichi,psichinormalized,n1,n2,n3) 2569 call copy_orbital(wan%psichi,psichinormalized,n1,n2,n3) 2570 2571 do isppol = 1,wan%nsppol 2572 do ispinor1 = 1,wan%nspinor 2573 do ispinor2 = 1,wan%nspinor 2574 index_l = 0 ! index_l is set to 0 at the beginning 2575 do iatom1 = 1,wan%natom_wan 2576 do il1 = 1,wan%nbl_atom_wan(iatom1) 2577 do im1 = 1,2*wan%latom_wan(iatom1)%lcalc(il1)+1 2578 index_l = index_l + 1 ! the line changes 2579 index_c = 1 ! counter_c is set to one each time the line changes 2580 do iatom2 = 1,wan%natom_wan 2581 do il2 = 1,wan%nbl_atom_wan(iatom2) 2582 do im2 = 1,2*wan%latom_wan(iatom2)%lcalc(il2)+1 2583 do ikpt = 1,wan%nkpt 2584 operwansquare(ikpt,isppol,index_l+wan%size_wan*(ispinor1-1),index_c+wan%size_wan*(ispinor2-1)) & 2585 & = operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2) 2586 end do 2587 index_c = index_c + 1 2588 end do !im2 2589 end do !il2 2590 end do ! iatom2 (the line changes) 2591 end do ! im1 2592 end do ! il1 2593 end do !iatom1 2594 end do 2595 end do 2596 end do 2597 2598 2599 2600 2601 !!Write the overlap matrix for ikpt = 1 in a nice shape 2602 ! do isppol = 1,wan%nsppol 2603 ! do il1 = 1,size(operwansquare,3) !! dummy variable without any meaning 2604 ! write(mat_writing,'(a,i0,i0)') 'Overlap matrix before orthonormalization 1 ',isppol,il1 2605 ! do il2 = 1,size(operwansquare,4) 2606 ! write(mat_writing2,'(F10.6)') real(operwansquare(1,isppol,il1,il2)) 2607 ! mat_writing = trim(mat_writing)//trim(mat_writing2) 2608 ! end do 2609 ! print*,trim(mat_writing) 2610 ! end do 2611 ! end do 2612 2613 2614 2615 !take the square root inverse of operwansquare for normalization purposes 2616 ABI_ALLOCATE(tmp_operwansquare,(wan%nspinor*wan%size_wan,wan%nspinor*wan%size_wan)) 2617 do isppol = 1,wan%nsppol 2618 do ikpt = 1,wan%nkpt 2619 tmp_operwansquare(:,:)=operwansquare(ikpt,isppol,:,:) 2620 call invsqrt_matrix(tmp_operwansquare,wan%nspinor*wan%size_wan) 2621 operwansquare(ikpt,isppol,:,:)=tmp_operwansquare(:,:) 2622 end do 2623 end do 2624 ABI_DEALLOCATE(tmp_operwansquare) 2625 2626 do ikpt = 1,wan%nkpt 2627 do iband = 1,wan%bandf_wan-wan%bandi_wan+1 2628 do iatom1 = 1,wan%natom_wan 2629 do il1 = 1,wan%nbl_atom_wan(iatom1) 2630 psichinormalized(ikpt,iband,iatom1)%atom(il1)%matl = czero 2631 end do 2632 end do 2633 end do 2634 end do 2635 2636 2637 2638 ! compute the new psichi normalized 2639 do isppol = 1,wan%nsppol 2640 do ispinor1 = 1,wan%nspinor 2641 do iband = 1,wan%bandf_wan-wan%bandi_wan+1 2642 index_l = 0 2643 do iatom1 = 1,wan%natom_wan 2644 do il1 = 1,wan%nbl_atom_wan(iatom1) 2645 do im1 = 1,2*wan%latom_wan(iatom1)%lcalc(il1)+1 2646 ! sum 2647 do ispinor2 = 1,wan%nspinor 2648 index_l = index_l + 1 ! the line changes 2649 index_c = 1 ! when the line changes, index_c is set to 1 2650 do iatom2 = 1,wan%natom_wan 2651 do il2 = 1,wan%nbl_atom_wan(iatom2) 2652 do im2 = 1,2*wan%latom_wan(iatom2)%lcalc(il2)+1 2653 do ikpt = 1,wan%nkpt 2654 psichinormalized(ikpt,iband,iatom1)%atom(il1)%matl(im1,isppol,ispinor1) =& 2655 psichinormalized(ikpt,iband,iatom1)%atom(il1)%matl(im1,isppol,ispinor1) +& 2656 wan%psichi(ikpt,iband,iatom2)%atom(il2)%matl(im2,isppol,ispinor2)*& 2657 operwansquare(ikpt,isppol,index_l+wan%size_wan*(ispinor1-1),index_c+& 2658 wan%size_wan*(ispinor2-1)) 2659 end do 2660 index_c = index_c + 1 2661 end do !im2 2662 end do !il2 2663 end do !iatom2 2664 end do ! ispinor2 2665 end do !im1 2666 end do ! il1 2667 end do ! iatom1 2668 end do ! iband 2669 end do ! ispinor1 2670 end do !isppol 2671 ! copy the new psichi normalized 2672 call copy_orbital(psichinormalized,wan%psichi,n1,n2,n3) 2673 call destroy_orbital(psichinormalized,n1,n2,n3) 2674 ABI_DATATYPE_DEALLOCATE(psichinormalized) 2675 2676 call destroy_operwan(wan,operwan) 2677 ABI_DATATYPE_DEALLOCATE(operwan) 2678 2679 2680 2681 2682 2683 2684 2685 2686 !! 2687 !! !------------------------------------------------------------- 2688 !! !check if the new norm is one 2689 2690 ABI_DATATYPE_ALLOCATE(operwan,(wan%nkpt,wan%natom_wan,wan%natom_wan)) 2691 call initialize_operwan(wan,operwan) 2692 do ikpt = 1,wan%nkpt 2693 call compute_oper_ks2wan(wan,operks,operwan,ikpt) 2694 end do 2695 2696 2697 do isppol = 1,wan%nsppol 2698 do ikpt = 1,wan%nkpt 2699 do iatom1 = 1,wan%natom_wan 2700 do iatom2 = 1,wan%natom_wan 2701 do il1 = 1,wan%nbl_atom_wan(iatom1) 2702 do il2 = 1,wan%nbl_atom_wan(iatom2) 2703 do im1 = 1,2*wan%latom_wan(iatom1)%lcalc(il1)+1 2704 do im2 = 1,2*wan%latom_wan(iatom2)%lcalc(il2)+1 2705 do ispinor1 = 1,wan%nspinor 2706 do ispinor2 = 1,wan%nspinor 2707 if (iatom1.eq.iatom2 .and. il1.eq.il2 .and. im1.eq.im2 .and. ispinor1.eq.ispinor2) then 2708 if (abs(cmplx(1.0,0.0,dpc)-& 2709 & operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2)) > 1d-8) then 2710 write(message,'(a,i0,a,F18.11)') 'Normalization error for ikpt =',ikpt,& 2711 & ' on diag, value = ',& 2712 & abs(operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2)) 2713 MSG_ERROR(message) 2714 end if 2715 else 2716 if (abs(operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2)) > 1d-8) then 2717 write(message,'(a,i0,a,F10.3)') 'Normalization error for ikpt =',ikpt,& 2718 & ' not on diag, value = ',& 2719 & abs(operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2)) 2720 MSG_ERROR(message) 2721 end if 2722 end if 2723 end do 2724 end do 2725 end do 2726 end do 2727 end do 2728 end do 2729 end do 2730 end do 2731 end do 2732 end do 2733 2734 2735 !!Uncomment to print the overlap matrix in the log file (for ikpt = 1) 2736 ! do isppol = 1,wan%nsppol 2737 ! do ispinor1 = 1,wan%nspinor 2738 ! do ispinor2 = 1,wan%nspinor 2739 ! index_l = 0 ! index_l is set to 0 at the beginning 2740 ! do iatom1 = 1,wan%natom_wan 2741 ! do il1 = 1,wan%nbl_atom_wan(iatom1) 2742 ! do im1 = 1,2*wan%latom_wan(iatom1)%lcalc(il1)+1 2743 ! index_l = index_l + 1 ! the line changes 2744 ! index_c = 1 ! counter_c is set to one each time the line changes 2745 ! do iatom2 = 1,wan%natom_wan 2746 ! do il2 = 1,wan%nbl_atom_wan(iatom2) 2747 ! do im2 = 1,2*wan%latom_wan(iatom2)%lcalc(il2)+1 2748 ! do ikpt = 1,wan%nkpt 2749 ! operwansquare(ikpt,isppol,index_l+wan%size_wan*(ispinor1-1),& 2750 ! &index_c+wan%size_wan*(ispinor2-1)) = operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,ispinor1,ispinor2) 2751 ! end do 2752 ! index_c = index_c + 1 2753 ! end do !im2 2754 ! end do !il2 2755 ! end do ! iatom2 (the line changes) 2756 ! end do ! im1 2757 ! end do ! il1 2758 ! end do !iatom1 2759 ! end do 2760 ! end do 2761 ! end do 2762 2763 !!Print the overlap matrix in a nice shape 2764 ! do isppol = 1,wan%nsppol 2765 ! do il1 = 1,size(operwansquare,3) !! dummy variable without any meaning 2766 ! write(mat_writing,'(a,i0,i0)') 'Overlap matrix after orthonormalization ',isppol,il1 2767 ! do il2 = 1,size(operwansquare,4) 2768 ! write(mat_writing2,'(F10.6)') real(operwansquare(9,isppol,il1,il2)) 2769 ! mat_writing = trim(mat_writing)//trim(mat_writing2) 2770 ! end do 2771 ! print*,trim(mat_writing) 2772 ! end do 2773 ! end do 2774 2775 2776 2777 2778 !! !---------------------------------------------------------------- 2779 ABI_DEALLOCATE(operwansquare) 2780 ABI_DEALLOCATE(operks) 2781 call destroy_operwan(wan,operwan) 2782 ABI_DATATYPE_DEALLOCATE(operwan) 2783 2784 2785 end subroutine normalization_plowannier
m_plowannier/operwan_realspace_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
operwan_realspace_type
FUNCTION
SOURCE
190 type, public :: operwan_realspace_type 191 192 type(operwan_type), allocatable :: position(:,:) 193 ! size (number of positions chosen for atom1, number of positions chosen for atom2) 194 195 end type operwan_realspace_type
m_plowannier/operwan_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
operwan_type
FUNCTION
SOURCE
173 type, public :: operwan_type 174 175 type(lorbital2_type), allocatable :: atom(:,:) 176 ! l chosen for each on of both atoms 177 178 end type operwan_type
m_plowannier/orbital_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
orbital_type
FUNCTION
SOURCE
133 type, public :: orbital_type 134 135 type(lorbital_type), allocatable :: atom(:) 136 ! details of the psichi coefficients for each atom 137 ! size of number of l chosen 138 139 end type orbital_type
m_plowannier/plowannier_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
plowannier_type
FUNCTION
SOURCE
207 type, public :: plowannier_type 208 209 integer :: nkpt 210 ! number of k points in Brillouin zone 211 212 integer :: bandi_wan 213 ! energy band minimum considered 214 215 integer :: bandf_wan 216 ! energy band maximum considered 217 218 integer :: natom_wan 219 ! number of atoms in a cell 220 221 integer :: size_wan 222 ! sum of all the m possible for every atom considered 223 224 integer, allocatable :: iatom_wan(:) 225 ! array of each atom's type 226 227 integer, allocatable :: nbl_atom_wan(:) 228 ! array of the number of l considered for each atom 229 230 type(latom_wan_type), allocatable :: latom_wan(:) 231 ! for each atom, it contains an array of the l we are interested in 232 233 integer, allocatable :: nbproj_atom_wan(:) 234 ! array of the number of projectors considered for each atom 235 236 type(projector_wan_type), allocatable :: projector_wan(:) 237 ! for each atom, it contains an array of the projectors we are interested in 238 239 type(position_wan_type), allocatable :: nposition(:) 240 ! array of the number of position considered for each atom 241 242 integer :: nsppol 243 ! number of polarization 244 245 integer :: nspinor 246 ! number of spinorial components 247 248 type(orbital_type), allocatable :: psichi(:,:,:) 249 ! arrays of psichi 250 251 integer, allocatable :: position(:,:) 252 ! size natom,3, gives the position of the cell for this atom (rprim coordinates) 253 254 real,allocatable :: kpt(:,:) 255 ! gives the coordinates in the BZ of the kpoint 256 ! size (3,nkpt) 257 258 real,allocatable :: wtk(:) 259 !weight of each kpoint 260 261 real,allocatable :: acell(:) 262 !size of the cell 263 264 end type plowannier_type
m_plowannier/position_wan_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
position_wan_type
FUNCTION
SOURCE
95 type, public :: position_wan_type 96 97 integer, allocatable :: pos(:,:) 98 ! size (number of position,3) 99 100 end type position_wan_type
m_plowannier/print_operwan [ Functions ]
[ Top ] [ m_plowannier ] [ Functions ]
NAME
print_operwan
FUNCTION
Print the Wannier operator (real space) in a latex file
INPUTS
wan, operwan, name
OUTPUT
PARENTS
m_plowannier
CHILDREN
SOURCE
2812 subroutine print_operwan(wan,operwan,name,convert) 2813 2814 !Arguments---------------------------------- 2815 2816 !This section has been created automatically by the script Abilint (TD). 2817 !Do not modify the following lines by hand. 2818 #undef ABI_FUNC 2819 #define ABI_FUNC 'print_operwan' 2820 !End of the abilint section 2821 2822 type(operwan_type),intent(in) :: operwan(:,:,:) 2823 type(plowannier_type), intent(in) :: wan 2824 character(len=*), intent(in) :: name 2825 real(dp), intent(in) :: convert 2826 2827 !Local variables---------------------------- 2828 integer :: iatom1,iatom2,pos1,pos2,il1,il2,im1,im2,isppol,ikpt,unt 2829 real(dp) :: sum 2830 character(len = 500) :: str1,str2,msg 2831 2832 if (open_file(name, msg, newunit=unt) /= 0) then 2833 MSG_ERROR(msg) 2834 end if 2835 2836 write(unt,'(a)') '\documentclass[11pt,a4paper,landscape]{article}' 2837 write(unt,'(a)') '\usepackage[T1]{fontenc}' 2838 write(unt,'(a)') '\usepackage{geometry,tabularx,graphicx}' 2839 write(unt,'(a)') '\geometry{left=0.5cm,right=0.5cm}' 2840 2841 write(unt,'(a)') '\begin{document}' 2842 write(unt,'(a)') '\noindent' 2843 2844 ! write(unt,'(a,i0,a,F7.3,a)') "% ",wan%natom_wan," atom in a ",wan%acell(1)," cell" 2845 ! write(unt,'(a)') "% atom isppol proj" 2846 2847 2848 do isppol = 1,wan%nsppol 2849 write(unt,'(a)') '\begin{figure}' 2850 write(unt,'(a)') '\resizebox{\linewidth}{!}{%' 2851 write(unt,'(a)') '$ \left(' 2852 2853 2854 write(str1,'(a)') '\begin{array}{' 2855 do iatom1 = 1,wan%natom_wan 2856 do pos1 = 1,size(wan%nposition(iatom1)%pos,1) 2857 do il1 = 1,wan%nbl_atom_wan(iatom1) 2858 do im1 = 1,2*wan%latom_wan(iatom1)%lcalc(il1)+1 2859 str1 = trim(str1)//"c" 2860 end do 2861 if (iatom1 .ne. wan%natom_wan .or. il1 .ne. wan%nbl_atom_wan(iatom1) .or. pos1 .ne. size(wan%nposition(iatom1)%pos,1)) then 2862 str1 = trim(str1)//'|' 2863 end if 2864 end do 2865 end do 2866 end do 2867 str1 = trim(str1)//'}' 2868 write(unt,'(a)') str1 2869 2870 2871 do iatom1 = 1,wan%natom_wan 2872 do pos1 = 1,size(wan%nposition(iatom1)%pos,1) 2873 do il1 = 1,wan%nbl_atom_wan(iatom1) 2874 do im1 = 1,2*wan%latom_wan(iatom1)%lcalc(il1)+1 2875 write(str1,'(a)') "" 2876 do iatom2 = 1,wan%natom_wan 2877 do pos2 = 1,size(wan%nposition(iatom2)%pos,1) 2878 do il2 = 1,wan%nbl_atom_wan(iatom2) 2879 do im2 = 1,2*wan%latom_wan(iatom2)%lcalc(il2)+1 2880 sum = 0 2881 do ikpt = 1,wan%nkpt 2882 sum = sum + convert*real(operwan(ikpt,iatom1,iatom2)%atom(il1,il2)%matl(im1,im2,isppol,1,1)*wan%wtk(ikpt)*& 2883 & exp(cmplx(0.0,1.0)*two_pi*(wan%kpt(1,ikpt)*( & 2884 & wan%nposition(iatom1)%pos(pos1,1)-wan%nposition(iatom2)%pos(pos2,1))+ & 2885 & wan%kpt(2,ikpt)*(wan%nposition(iatom1)%pos(pos1,2)-wan%nposition(iatom2)%pos(pos2,2))+ & 2886 & wan%kpt(3,ikpt)*(wan%nposition(iatom1)%pos(pos1,3)-wan%nposition(iatom2)%pos(pos2,3))))) 2887 end do 2888 write(str2,'(F10.6)') real(sum) 2889 if ( len_trim(str1) .ge. 2) then 2890 str1 = trim(str1)//"&"//trim(str2) 2891 else 2892 str1 = trim(str2) 2893 end if 2894 if (iatom2 .eq. wan%natom_wan .and. il2 .eq. wan%nbl_atom_wan(iatom2) .and. im2 & 2895 & .eq. 2*wan%latom_wan(iatom2)%lcalc(il2)+1 .and. pos2 .eq. size(wan%nposition(iatom2)%pos,1)) then 2896 str1 = trim(str1)//'\\' 2897 end if 2898 end do 2899 end do 2900 end do 2901 end do 2902 write(unt,'(a)') trim(str1) 2903 end do 2904 if (iatom1 .ne. wan%natom_wan .or. il1 .ne. wan%nbl_atom_wan(iatom1) .or. pos1 .ne. size(wan%nposition(iatom1)%pos,1)) then 2905 write(unt,'(a)') '\hline' 2906 end if 2907 end do 2908 end do 2909 end do 2910 2911 2912 2913 write(unt,'(a)') '\end{array} \right) $ }' 2914 2915 if (name(len_trim(name)-2:len(trim(name))) == 'gen') then 2916 write(unt,'(a,i0,a,F7.3,a)') "\caption{Energy matrix in real space for isppol = ", & 2917 & isppol," in a ",wan%acell(1), " a.u. cell}" 2918 end if 2919 2920 if (name(len_trim(name)-2:len(trim(name))) == 'occ') then 2921 write(unt,'(a,i0,a,F7.3,a)') "\caption{Occupation matrix in real space for isppol = ",& 2922 & isppol," in a ",wan%acell(1), " a.u. cell}" 2923 end if 2924 2925 2926 2927 write(unt,'(a)') '\end{figure}' 2928 write(unt,'(a)') '\end{document}' 2929 2930 end do 2931 close(unt) 2932 2933 2934 2935 end subroutine print_operwan 2936 2937 END MODULE m_plowannier
m_plowannier/projector_wan_type [ Types ]
[ Top ] [ m_plowannier ] [ Types ]
NAME
projector_wan_type
FUNCTION
SOURCE
77 type, public :: projector_wan_type 78 79 integer, allocatable :: lproj(:) 80 ! gives the list of the projector chosen 81 82 end type projector_wan_type