TABLE OF CONTENTS
- ABINIT/m_bseinterp
- m_bseinterp/int_alloc_work
- m_bseinterp/int_compute_overlaps
- m_bseinterp/int_free_work
- m_bseinterp/int_preprocess_tables
- m_bseinterp/interpolator_free
- m_bseinterp/interpolator_init
- m_bseinterp/interpolator_normalize
- m_haydock/int_compute_corresp
- m_haydock/interpolator_t
ABINIT/m_bseinterp [ Modules ]
NAME
m_bseinterp
FUNCTION
COPYRIGHT
Copyright (C) 2014-2024 ABINIT group (M.Giantomassi, Y. Gillet) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
15 #if defined HAVE_CONFIG_H 16 #include "config.h" 17 #endif 18 19 #include "abi_common.h" 20 21 MODULE m_bseinterp 22 23 use defs_basis 24 use m_abicore 25 use m_bs_defs 26 use m_xmpi 27 use m_errors 28 use m_nctk 29 use m_haydock_io 30 use m_linalg_interfaces 31 use netcdf 32 33 use m_fstrings, only : indent, strcat, sjoin, itoa 34 use defs_datatypes, only : pseudopotential_type 35 use m_hide_blas, only : xdotc 36 use m_fft_mesh, only : calc_ceigr 37 use m_crystal, only : crystal_t 38 use m_bz_mesh, only : kmesh_t 39 use m_double_grid, only : double_grid_t, get_kpt_from_indices_coarse 40 use m_wfd, only : wfdgw_t 41 use m_pawtab, only : pawtab_type 42 43 implicit none 44 45 private
m_bseinterp/int_alloc_work [ Functions ]
[ Top ] [ m_bseinterp ] [ Functions ]
NAME
int_alloc_work
FUNCTION
Allocate temporary arrays
INPUTS
OUTPUT
SOURCE
216 subroutine int_alloc_work(interpolator, work_size) 217 218 !Arguments --------------------------- 219 !scalars 220 integer,intent(in) :: work_size 221 type(interpolator_t),intent(inout) :: interpolator 222 223 !***************************************************************************** 224 225 ABI_MALLOC(interpolator%btemp,(work_size)) 226 ABI_MALLOC(interpolator%ctemp,(work_size)) 227 228 end subroutine int_alloc_work
m_bseinterp/int_compute_overlaps [ Functions ]
[ Top ] [ m_bseinterp ] [ Functions ]
NAME
int_compute_overlaps
FUNCTION
Compute the overlaps prefactors
INPUTS
OUTPUT
SOURCE
275 subroutine int_compute_overlaps(interpolator, double_grid, Wfd_dense, Wfd_coarse, & 276 & Kmesh_dense, Kmesh_coarse, BSp, Cryst, Psps, Pawtab) 277 278 !Arguments --------------------------- 279 !scalars 280 type(interpolator_t),intent(inout) :: interpolator 281 type(double_grid_t),intent(in),target :: double_grid 282 type(wfdgw_t),intent(inout) :: Wfd_dense, Wfd_coarse 283 type(kmesh_t),intent(in) :: Kmesh_dense, Kmesh_coarse 284 type(excparam),intent(in) :: BSp 285 type(crystal_t),intent(in) :: Cryst 286 type(pseudopotential_type),intent(in) :: Psps 287 !arrays 288 type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd_coarse%usepaw) 289 290 !Local variables --------------------- 291 !scalars 292 integer :: nprocs, my_rank 293 integer :: ierr 294 integer :: nfft, nspinor, nsppol, nvert 295 integer :: ib_coarse, ib_dense 296 integer :: ik_coarse, ik_dense 297 integer :: spin 298 integer :: iorder 299 integer :: ivertex, ix, iy, iz 300 integer :: bstart, bstop 301 real(dp),parameter :: threshold = 0.1_dp 302 complex(gwpc) :: ovlp 303 !arrays 304 integer :: curindices_dense(6), curindices_coarse(3) 305 integer :: neighbour(3) 306 integer :: g0(3),g01(3),diffg0(3) 307 complex(gwpc),allocatable :: ur_coarse(:),ur_dense(:) 308 complex(gwpc),allocatable :: ceigr(:) 309 !arrays 310 311 !***************************************************************************** 312 313 nprocs = xmpi_comm_size(Wfd_coarse%comm) 314 my_rank = xmpi_comm_rank(Wfd_coarse%comm) 315 316 ABI_UNUSED(Pawtab(1)%basis_size) 317 318 ! Ensure Wfd and Wfd_coarse use the same FFT mesh. 319 call wfd_dense%change_ngfft(Cryst,Psps,Wfd_coarse%ngfft) 320 nfft = Wfd_coarse%nfft 321 nspinor = Wfd_coarse%nspinor 322 nsppol = Bsp%nsppol 323 nvert = interpolator%nvert 324 325 ! Allocate workspace for wavefunctions in real space. 326 ABI_MALLOC(ur_coarse,(nfft*nspinor)) 327 ABI_MALLOC(ur_dense,(nfft*nspinor)) 328 ABI_MALLOC(ceigr,(nfft*nspinor)) 329 330 interpolator%overlaps = czero 331 332 ! TODO 333 ! 1) Choose whether we want to compute only dvv, dcc or all dbb 334 ! 2) Check the ordering of the loops 335 ! 3) Improve vertex -> neighbour (in double_grid ?) 336 do spin = 1,nsppol 337 do ik_dense = 1,double_grid%nbz_dense 338 339 ! MPI parallelization 340 ! We assume that each node owns in memory the full set of wavefunctions 341 ! both coarse and dense k-mesh and both spins. 342 if (mod(ik_dense, nprocs) /= my_rank) cycle 343 344 ! From ik_dense -> indices_dense 345 iorder = double_grid%iktoint_dense(ik_dense) 346 g01 = double_grid%g0_dense(:,iorder) 347 curindices_dense = double_grid%indices_dense(:,iorder) 348 349 do ivertex = 1,nvert 350 351 ! From vertex to neighbour 352 ! TODO improve this part + permit to choose other neighbour (e.g. nearest neighbour for RL) 353 if(nvert > 1) then 354 ix = (ivertex-1)/4 355 iy = (ivertex-ix*4-1)/2 356 iz = (ivertex-ix*4-iy*2-1) 357 else 358 ix = (BSp%rl_nb-1)/4 359 iy = (BSp%rl_nb-ix*4-1)/2 360 iz = (BSp%rl_nb-ix*4-iy*2-1) 361 end if 362 363 neighbour = [ix,iy,iz] 364 365 ! From indices_dense -> indices_coarse 366 curindices_coarse = curindices_dense(1:3) + neighbour(:) 367 368 ! From indices_coarse -> ik_ibz in the coarse mesh 369 call get_kpt_from_indices_coarse(curindices_coarse,double_grid%maxcomp_coarse,& 370 & double_grid%inttoik_coarse,double_grid%g0_coarse,double_grid%nbz_closedcoarse,ik_coarse,g0) 371 372 ! Take into account a possible umklapp between k_dense and k_coarse 373 diffg0 = g0 - g01 374 375 if (ANY(diffg0/=0)) then 376 ! WARNING works only with nspinor = 1 !!! 377 call calc_ceigr(diffg0,nfft,nspinor,Wfd_coarse%ngfft,ceigr) 378 end if 379 380 do ib_dense = BSp%lomo_spin(spin), BSp%humo_spin(spin) 381 ! ur(ib_dense, ik_dense) 382 call wfd_dense%sym_ur(Cryst,Kmesh_dense,ib_dense,ik_dense,spin,ur_dense) 383 384 if (ANY(diffg0/=0)) then 385 !ur_kbz = ur_kbz*e(ig0r) 386 ur_dense(:) = ur_dense(:)*ceigr(:) 387 end if 388 389 ! Uncomment for the complete overlap 390 !bstart = BSp%lomo_spin(spin); bstop = BSp%humo_spin(spin) 391 392 ! Compute only dvv or dcc 393 if (ib_dense <= BSp%homo_spin(spin)) then 394 ! if ib_dense is a valence band => loop on valence bands 395 bstart = BSp%lomo_spin(spin); bstop = BSp%homo_spin(spin) 396 else 397 ! if ib_dense is a conduction band => loop on conduction bands 398 bstart = BSp%lumo_spin(spin); bstop = BSp%humo_spin(spin) 399 end if 400 401 do ib_coarse = bstart, bstop 402 ! ur(ib_coarse, ik_coarse) 403 call wfd_coarse%sym_ur(Cryst,Kmesh_coarse,ib_coarse,ik_coarse,spin,ur_coarse) 404 405 ! ovlp = < u_{ib_coarse,ik_coarse} | u_{ib_dense,ik_dense} > 406 ovlp = xdotc(nfft,ur_coarse,1,ur_dense,1)/nfft 407 408 ! Filter too low values 409 if (ABS(ovlp) < threshold) ovlp = czero 410 411 interpolator%overlaps(ib_coarse,ib_dense,ivertex,ik_dense,spin) = ovlp 412 end do ! ib_coarse 413 414 !DBYG 415 ! write(std_out,*) "nb = ",neighbour 416 ! write(std_out,*) "(i1,i2,i3,j1,j2,j3) = ",curindices_dense 417 ! write(std_out,*) "ib = ",ib_dense 418 ! write(std_out,*) "Sum of dbb = ",REAL(SUM(GWPC_CONJG(interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin))*interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin))); call flush(std_out) 419 !ENDDBYG 420 421 end do ! ib_dense 422 end do ! ivertex 423 end do ! ik_dense 424 end do ! spin 425 426 ABI_FREE(ur_coarse) 427 ABI_FREE(ur_dense) 428 ABI_FREE(ceigr) 429 430 ! Gather results on each node. 431 call xmpi_sum(interpolator%overlaps,Wfd_coarse%comm,ierr) 432 433 end subroutine int_compute_overlaps
m_bseinterp/int_free_work [ Functions ]
[ Top ] [ m_bseinterp ] [ Functions ]
NAME
int_free_work
FUNCTION
Deallocate temporary arrays
INPUTS
OUTPUT
SOURCE
246 subroutine int_free_work(interpolator) 247 248 !Arguments --------------------------- 249 !scalars 250 type(interpolator_t),intent(inout) :: interpolator 251 252 !***************************************************************************** 253 254 ABI_SFREE(interpolator%btemp) 255 ABI_SFREE(interpolator%ctemp) 256 257 end subroutine int_free_work
m_bseinterp/int_preprocess_tables [ Functions ]
[ Top ] [ m_bseinterp ] [ Functions ]
NAME
int_preprocess_tables
FUNCTION
Pre-process tables to improve interpolation technique
INPUTS
OUTPUT
SOURCE
451 subroutine int_preprocess_tables(interpolator,double_grid) 452 453 !Argument ------------------------------------ 454 !scalars 455 type(interpolator_t),intent(inout) :: interpolator 456 type(double_grid_t),intent(in) :: double_grid 457 !arrays 458 459 !Local variables ----------------------------- 460 !scalars 461 integer :: iorder,ik_dense,ik_coarse 462 integer :: ix,iy,iz,ineighbour,curdim, curj 463 real(dp) :: interp_factor 464 !arrays 465 integer :: allxyz(3),curindices_dense(6) 466 integer,allocatable :: curindex(:) 467 468 !********************************************* 469 ABI_MALLOC(curindex,(double_grid%nbz_coarse)) 470 curindex = 1 471 472 interpolator%interp_factors = zero 473 474 do ik_dense = 1,double_grid%nbz_dense 475 476 ! From ik_ibz in the dense mesh -> indices_dense 477 iorder = double_grid%iktoint_dense(ik_dense) 478 !g01 = double_grid%g0_dense(:,iorder) 479 480 ! From indices_dense -> indices_coarse 481 curindices_dense = double_grid%indices_dense(:,iorder) 482 483 ik_coarse = double_grid%dense_to_coarse(ik_dense) 484 485 ! Compute multi-linear interpolation factors 486 ! Loop over the neighbours 487 do ineighbour = 1,interpolator%nvert 488 !TODO helper function from [ix,iy,iz] -> ineighbour and vice versa 489 ix = (ineighbour-1)/4 490 iy = (ineighbour-ix*4-1)/2 491 iz = (ineighbour-ix*4-iy*2-1) 492 allxyz = [ix,iy,iz] 493 interp_factor = one 494 do curdim = 1,3 495 if (interpolator%method == BSE_INTERP_RL) then 496 cycle 497 else if(interpolator%method == BSE_INTERP_RL2) then 498 if (curdim /= 3) cycle 499 end if 500 curj = curindices_dense(3+curdim) 501 interp_factor = interp_factor*((allxyz(curdim)*(curj*1.0/double_grid%kmult(curdim)))& 502 & +((1-allxyz(curdim))*(1-(curj*1.0/double_grid%kmult(curdim))))) 503 end do 504 interpolator%interp_factors(ineighbour,curindex(ik_coarse)) = interp_factor 505 end do 506 507 curindex(ik_coarse) = curindex(ik_coarse) + 1 508 end do 509 510 ABI_FREE(curindex) 511 512 end subroutine int_preprocess_tables
m_bseinterp/interpolator_free [ Functions ]
[ Top ] [ m_bseinterp ] [ Functions ]
NAME
interpolator_free
FUNCTION
Destroy the interpolator object in memory
INPUTS
OUTPUT
SOURCE
674 subroutine interpolator_free(interpolator) 675 676 !Arguments --------------------------- 677 type(interpolator_t),intent(inout) :: interpolator 678 679 !***************************************************************************** 680 681 ABI_SFREE(interpolator%overlaps) 682 ABI_SFREE(interpolator%corresp) 683 ABI_SFREE(interpolator%interp_factors) 684 if( associated(interpolator%double_grid) ) nullify(interpolator%double_grid) 685 686 end subroutine interpolator_free
m_bseinterp/interpolator_init [ Functions ]
[ Top ] [ m_bseinterp ] [ Functions ]
NAME
interpolator_init
FUNCTION
Construct the interpolator object
INPUTS
OUTPUT
SOURCE
123 subroutine interpolator_init(interpolator, double_grid, Wfd_dense, Wfd_coarse, & 124 & Kmesh_dense, Kmesh_coarse, BSp, Cryst, Psps, Pawtab, method) 125 126 !Arguments --------------------------- 127 !scalars 128 integer,intent(in) :: method 129 type(interpolator_t),intent(inout) :: interpolator 130 type(double_grid_t),intent(in),target :: double_grid 131 type(wfdgw_t),intent(inout) :: Wfd_dense, Wfd_coarse 132 type(kmesh_t),intent(in) :: Kmesh_dense, Kmesh_coarse 133 type(excparam),intent(in) :: BSp 134 type(crystal_t),intent(in) :: Cryst 135 type(pseudopotential_type),intent(in) :: Psps 136 !arrays 137 type(pawtab_type),intent(in) :: Pawtab(Cryst%ntypat*Wfd_coarse%usepaw) 138 139 !Local variables --------------------- 140 !scalars 141 integer :: nsppol, nvert 142 integer :: maxnreh, nreh1, nreh2 143 integer :: mbandc, mbandd, nbzd 144 real(dp),parameter :: threshold = 0.1_dp 145 !arrays 146 character(len=500) :: msg 147 148 !***************************************************************************** 149 150 ABI_CHECK(Wfd_coarse%usepaw==0, "PAW not yet supported") 151 ABI_CHECK(BSp%nsppol==1, "nsppol != 1 not yet implemented") 152 ABI_CHECK(Wfd_coarse%nspinor==1, "nspinor != 1 not supported") 153 154 ABI_UNUSED(Pawtab(1)%basis_size) 155 !paw_overlap(cprj1,cprj2,typat,pawtab,spinor_comm) result(onsite) 156 157 interpolator%double_grid => double_grid 158 interpolator%mband_dense = Wfd_dense%mband 159 interpolator%mband_coarse = Wfd_coarse%mband 160 interpolator%method = method 161 interpolator%nsppol = BSp%nsppol 162 163 SELECT CASE(method) 164 CASE (BSE_INTERP_YG) 165 nvert = 8 166 CASE (BSE_INTERP_RL2) 167 nvert = 2 168 CASE (BSE_INTERP_RL) 169 nvert = 1 170 CASE DEFAULT 171 write(msg,'(a,i0)') "Wrong interpolation method: ",method 172 ABI_ERROR(msg) 173 END SELECT 174 175 interpolator%nvert = nvert 176 177 mbandc = interpolator%mband_coarse 178 mbandd = interpolator%mband_dense 179 nbzd = double_grid%nbz_dense 180 nsppol = interpolator%nsppol 181 ABI_MALLOC(interpolator%overlaps,(mbandc,mbandd,nvert,nbzd,nsppol)) 182 183 call int_compute_overlaps(interpolator,double_grid, Wfd_dense, Wfd_coarse, Kmesh_dense, & 184 & Kmesh_coarse, BSp, Cryst, Psps, Pawtab) 185 186 ABI_MALLOC(interpolator%interp_factors,(nvert,double_grid%ndiv)) 187 188 call int_preprocess_tables(interpolator,double_grid) 189 190 nreh1 = BSp%nreh(1) 191 nreh2 = nreh1; if(BSp%nsppol == 2) nreh2 = BSp%nreh(2) 192 maxnreh = MAX(nreh1,nreh2) 193 194 ABI_MALLOC(interpolator%corresp,(maxnreh,interpolator%nvert,interpolator%nsppol)) 195 196 call int_compute_corresp(interpolator,BSp,double_grid) 197 198 end subroutine interpolator_init
m_bseinterp/interpolator_normalize [ Functions ]
[ Top ] [ m_bseinterp ] [ Functions ]
NAME
interpolator_normalize
FUNCTION
Normalize the overlaps so that \sum_{ib} | d_{kk'}^{b,ib} | ^2 = 1
INPUTS
OUTPUT
SOURCE
620 subroutine interpolator_normalize(interpolator) 621 622 !Arguments --------------------------- 623 !scalars 624 type(interpolator_t),intent(inout) :: interpolator 625 626 !Local variables --------------------- 627 !scalars 628 integer :: spin, ivertex, ib_dense, ik_dense 629 complex(gwpc) :: sum_ovlp 630 !arrays 631 complex(gwpc),allocatable :: overlaps(:) 632 633 !***************************************************************************** 634 635 ABI_MALLOC(overlaps,(interpolator%mband_coarse)) 636 do spin = 1, interpolator%nsppol 637 do ivertex = 1, interpolator%nvert 638 do ib_dense = 1, interpolator%mband_dense 639 do ik_dense = 1, interpolator%double_grid%nbz_dense 640 overlaps(:) = interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin) 641 sum_ovlp = SQRT(REAL(SUM(GWPC_CONJG(overlaps(:))*overlaps(:)))) 642 if (ABS(sum_ovlp) > tol6) then 643 overlaps(:) = overlaps(:)/sum_ovlp 644 else 645 overlaps(:) = czero 646 end if 647 interpolator%overlaps(:,ib_dense,ivertex,ik_dense,spin) = overlaps(:) 648 end do 649 end do 650 end do 651 end do 652 ABI_FREE(overlaps) 653 654 end subroutine interpolator_normalize
m_haydock/int_compute_corresp [ Functions ]
[ Top ] [ m_haydock ] [ Functions ]
NAME
int_compute_corresp
FUNCTION
INPUTS
BSp<type(excparam)=The parameter for the Bethe-Salpeter run. grid <double_grid_t>=Correspondence between coarse and fine k-grid spin=Spin index.
OUTPUT
corresp(Bsp%nreh(spin),8)= Correspondence between a transition on the coarse mesh and its i-th neighbour for i in [1,2,..,8]. TODO: Some operations are faster if we allocate with shape (8,nreh(spin))
SOURCE
537 subroutine int_compute_corresp(interpolator,BSp,double_grid) 538 539 !Arguments ------------------------------------ 540 !scalars 541 type(excparam),intent(in) :: BSp 542 type(double_grid_t),intent(in) :: double_grid 543 type(interpolator_t),intent(inout) :: interpolator 544 545 !Local variables ------------------------------ 546 !scalars 547 integer :: spin 548 integer :: itt,ik_dense,ik_coarse,iorder,it_coarse 549 integer :: ic,iv,ik_coarse0,it_coarse0,iovlp,ix,iy,iz 550 !arrays 551 integer :: curindices_dense(6),curindices_coarse(3),g0(3),g01(3),neighbour(3) 552 553 !************************************************************************ 554 555 do spin=1,interpolator%nsppol 556 do itt=1,BSp%nreh_interp(spin) 557 ! From dense itt -> ik_dense, ic, iv 558 ik_dense = BSp%Trans_interp(itt,spin)%k 559 ic = BSp%Trans_interp(itt,spin)%c 560 iv = BSp%Trans_interp(itt,spin)%v 561 562 ! From ik_dense -> indices_dense 563 iorder = double_grid%iktoint_dense(ik_dense) 564 g01 = double_grid%g0_dense(:,iorder) 565 566 ! Index of the k-point in the coarse mesh. 567 ik_coarse0 = double_grid%dense_to_coarse(ik_dense) 568 it_coarse0 = BSp%vcks2t(iv,ic,ik_coarse0,spin) 569 570 ! From indices_dense -> indices_coarse 571 curindices_dense = double_grid%indices_dense(:,iorder) 572 573 ! Loop over the 8 neighbors. 574 do iovlp = 1,interpolator%nvert 575 576 !TODO : helper function from [ix,iy,iz] -> iovlp and vice versa 577 if(interpolator%nvert > 1) then 578 ix = (iovlp-1)/4 579 iy = (iovlp-ix*4-1)/2 580 iz = (iovlp-ix*4-iy*2-1) 581 else 582 ix = (BSp%rl_nb-1)/4 583 iy = (BSp%rl_nb-ix*4-1)/2 584 iz = (BSp%rl_nb-ix*4-iy*2-1) 585 end if 586 neighbour = [ix,iy,iz] 587 588 curindices_coarse = curindices_dense(1:3) + neighbour(:) 589 590 ! From indices_coarse -> ik_ibz in the coarse mesh 591 call get_kpt_from_indices_coarse(curindices_coarse,double_grid%maxcomp_coarse,& 592 & double_grid%inttoik_coarse,double_grid%g0_coarse,double_grid%nbz_closedcoarse,ik_coarse,g0) 593 594 ! From ik_coarse, ic, iv to it_coarse 595 it_coarse = BSp%vcks2t(iv,ic,ik_coarse,spin) 596 597 interpolator%corresp(it_coarse0,iovlp,spin) = it_coarse 598 end do 599 end do ! itt 600 end do 601 602 end subroutine int_compute_corresp
m_haydock/interpolator_t [ Types ]
[ Top ] [ m_haydock ] [ Types ]
NAME
interpolator_t
FUNCTION
Store the overlap matrix elements needed for the interpolation of the BSE Hamiltonian
TODO
Decide if we want to make the number of bands k-dependent.
SOURCE
62 type,public :: interpolator_t 63 64 integer :: nvert=8 65 ! Number of vertices for interpolation 66 67 integer :: method 68 ! Interpolation method (YG or Rohlfing & Louie or ...) 69 70 integer :: mband_dense, mband_coarse 71 ! Max number of bands dense and coarse 72 73 integer :: nsppol 74 ! Number of spin channels 75 76 integer, allocatable :: corresp(:,:,:) 77 ! corresp(max_nreh,nvert,spin) 78 ! it_coarse, idiv -> it_coarse (idiv-th neighbour) 79 80 real(dp),allocatable :: interp_factors(:,:) 81 ! interp_factors(nvert,ndiv) 82 ! index_in_fine_box -> k-point in Trans_interp 83 84 complex(gwpc),allocatable :: overlaps(:,:,:,:,:) 85 ! Overlaps between dense and coarse mesh 86 ! overlaps(mband_coarse,mband_dense,ivertex_coarse,double_grid%nkpt_dense,spin) 87 88 complex(dpc),allocatable :: btemp(:), ctemp(:) 89 ! Temporary arrays for work 90 91 ! Pointers to datatypes that are already in memory 92 type(double_grid_t),pointer :: double_grid => null() 93 ! Mapping between coarse and dense mesh 94 95 end type interpolator_t 96 97 public :: interpolator_init ! Construct the object 98 public :: interpolator_free ! Free memory 99 public :: interpolator_normalize ! Normalize the overlaps 100 public :: int_alloc_work ! Alloc temp memory 101 public :: int_free_work ! Free temp memory