TABLE OF CONTENTS
- ABINIT/gspot_transgrid_and_pack
- ABINIT/m_hamiltonian
- m_hamiltonian/copy_hamiltonian
- m_hamiltonian/destroy_hamiltonian
- m_hamiltonian/destroy_rf_hamiltonian
- m_hamiltonian/gs_hamiltonian_type
- m_hamiltonian/init_hamiltonian
- m_hamiltonian/init_rf_hamiltonian
- m_hamiltonian/load_k_hamiltonian
- m_hamiltonian/load_k_rf_hamiltonian
- m_hamiltonian/load_kprime_hamiltonian
- m_hamiltonian/load_spin_hamiltonian
- m_hamiltonian/load_spin_rf_hamiltonian
- m_hamiltonian/pawdij2e1kb
- m_hamiltonian/pawdij2ekb
- m_hamiltonian/rf_hamiltonian_type
ABINIT/gspot_transgrid_and_pack [ Functions ]
NAME
gspot_transgrid_and_pack
FUNCTION
Set up local potential vlocal on the coarse FFT mesh with proper dimensioning from vtrial given on the fine mesh. Also take into account the spin.
INPUTS
isppol=Spin polarization. usepaw=1 if PAW paral_kgb: 1 if paral_kgb nfft=(effective) number of FFT grid points on the coarse mesh (for this processor) ngfft(18)contain all needed information about 3D FFT, for the coarse FFT mesh. see ~abinit/doc/variables/vargs.htm#ngfft nfftf=(effective) number of FFT grid points on the fine mesh (for this processor) nvloc==1 if nspden <=2, nvloc==4 for nspden==4, nspden=Number of spin density components. ncomp=Number of extra components in vtrial and vlocal (e.g. 1 if LDA/GGA pot, 4 for Meta-GGA, etc). pawfgr <type(pawfgr_type)>=fine grid parameters and related data vtrial(nfftf,nspden)=INPUT potential Vtrial(r). mpi_enreg=information about MPI parallelization
OUTPUT
vlocal(n4,n5,n6,nvloc,ncomp): Potential on the coarse grid.
SIDE EFFECTS
SOURCE
2062 subroutine gspot_transgrid_and_pack(isppol, usepaw, paral_kgb, nfft, ngfft, nfftf, & 2063 nspden, nvloc, ncomp, pawfgr, mpi_enreg, vtrial, vlocal) 2064 2065 !Arguments ------------------------------- 2066 integer,intent(in) :: isppol, nspden, ncomp, usepaw, paral_kgb, nfft, nfftf, nvloc 2067 type(pawfgr_type), intent(in) :: pawfgr 2068 type(MPI_type), intent(in) :: mpi_enreg 2069 !arrays 2070 integer,intent(in) :: ngfft(18) 2071 real(dp),intent(inout) :: vtrial(nfftf, nspden, ncomp) 2072 real(dp),intent(out) :: vlocal(ngfft(4), ngfft(5), ngfft(6), nvloc, ncomp) 2073 2074 2075 !Local variables------------------------------- 2076 !scalars 2077 integer :: n1,n2,n3,n4,n5,n6,ispden,ic 2078 real(dp) :: rhodum(1) 2079 real(dp),allocatable :: cgrvtrial(:,:), vlocal_tmp(:,:,:) 2080 2081 ! ************************************************************************* 2082 2083 ! Coarse mesh. 2084 n1=ngfft(1); n2=ngfft(2); n3=ngfft(3) 2085 n4=ngfft(4); n5=ngfft(5); n6=ngfft(6) 2086 2087 ! Set up local potential vlocal with proper dimensioning, from vtrial 2088 ! Also take into account the spin. 2089 if (nspden /= 4) then 2090 if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then 2091 ! Fine mesh == Coarse mesh. 2092 do ic=1,ncomp 2093 call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal(:,:,:,:,ic),2) 2094 end do 2095 else 2096 ! Transfer from fine mesh to coarse and then pack data 2097 ABI_MALLOC(cgrvtrial,(nfft,nspden)) 2098 do ic=1,ncomp 2099 call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,& 2100 rhodum,rhodum,cgrvtrial,vtrial(:,:,ic)) 2101 call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,& 2102 cgrvtrial,vlocal(:,:,:,:,ic),2) 2103 end do 2104 ABI_FREE(cgrvtrial) 2105 end if 2106 else 2107 ! nspden == 4. replace isppol by loop over ispden. 2108 ABI_MALLOC(vlocal_tmp, (n4,n5,n6)) 2109 if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then 2110 ! Fine mesh == Coarse mesh. 2111 do ic=1,ncomp 2112 do ispden=1,nspden 2113 call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal_tmp,2) 2114 vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:) 2115 end do 2116 end do 2117 else 2118 ! Transfer from fine mesh to coarse and then pack data 2119 ABI_MALLOC(cgrvtrial,(nfft,nspden)) 2120 do ic=1,ncomp 2121 call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,rhodum,rhodum,cgrvtrial,vtrial(:,:,ic)) 2122 do ispden=1,nspden 2123 call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,cgrvtrial,vlocal_tmp,2) 2124 vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:) 2125 end do 2126 end do 2127 ABI_FREE(cgrvtrial) 2128 end if 2129 ABI_FREE(vlocal_tmp) 2130 end if ! nspden 2131 2132 end subroutine gspot_transgrid_and_pack
ABINIT/m_hamiltonian [ Modules ]
NAME
m_hamiltonian
FUNCTION
This module provides the definition of the gs_hamiltonian_type and rf_hamiltonian_type datastructures used in the "getghc" and "getgh1c" routines to apply the Hamiltonian (or its derivative) on a wavefunction. Methods to initialize or destroy the objects are defined here.
TODO
All array pointers in H datatypes should be declared as contiguous for efficient reasons (well, here performance is critical) Client code should make sure they always point contiguous targets.
COPYRIGHT
Copyright (C) 2009-2024 ABINIT group (MG, MT) 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
24 #if defined HAVE_CONFIG_H 25 #include "config.h" 26 #endif 27 28 #include "abi_common.h" 29 30 module m_hamiltonian 31 32 use iso_fortran_env, only : int32,int64,real32,real64 33 34 use defs_basis 35 use m_abicore 36 use m_errors 37 use m_xmpi 38 39 use defs_datatypes, only : pseudopotential_type 40 use defs_abitypes, only : MPI_type 41 use m_copy, only : addr_copy 42 use m_geometry, only : metric 43 use m_pawtab, only : pawtab_type 44 use m_pawfgr, only : pawfgr_type 45 use m_fftcore, only : sphereboundary 46 use m_fft, only : fftpac 47 use m_fourier_interpol, only : transgrid 48 use m_pawcprj, only : pawcprj_getdim 49 use m_paw_ij, only : paw_ij_type 50 use m_paral_atom, only : get_my_atmtab, free_my_atmtab 51 use m_electronpositron, only : electronpositron_type, electronpositron_calctype 52 use m_kg, only : ph1d3d, getph 53 use m_fock, only : fock_common_type, fock_BZ_type, fock_ACE_type, fock_type 54 55 #if defined HAVE_GPU_CUDA 56 use m_manage_cuda 57 #endif 58 59 #if defined HAVE_FC_ISO_C_BINDING 60 use, intrinsic :: iso_c_binding, only : c_ptr,c_loc,c_f_pointer,c_int32_t,c_int64_t,c_size_t 61 #endif 62 63 #if defined HAVE_GPU && defined HAVE_YAKL 64 use gator_mod 65 #endif 66 67 implicit none 68 69 private 70 71 public :: pawdij2ekb 72 public :: pawdij2e1kb 73 public :: gspot_transgrid_and_pack ! Set up local potential vlocal on the coarse FFT mesh from vtrial on the fine mesh. 74 75 ! These constants select how H is applied in reciprocal space 76 integer,parameter,public :: KPRIME_H_K=1, K_H_KPRIME=2, K_H_K=3, KPRIME_H_KPRIME=4
m_hamiltonian/copy_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
copy_hamiltonian
INPUTS
gs_hamk_in<gs_hamiltonian_type>=Structured datatype completely initialized, to be copied.
FUNCTION
Copy a gs_hamiltonian_type variable (gs_hamk_in) in another (gs_hamk_out). In contrast to an assignment statement (gs_hamk_out=gs_hamk_in), this subroutine allocate memory space for the pointers contained in the data structure (gs_hamk_out) and copy the content of the corresponding memory space of gs_hamk_in in it. In contrast, the assignment statement would only associate the pointers of gs_hamk_out to the same memory space than the corresponding ones in gs_hamk_in. This can cause trouble if one data structure is destroyed before a reading/writing statement for the other structure, causing access to unallocated memory space (silently, without segmentation fault being generated).
OUTPUT
gs_hamk_out<gs_hamiltonian_type>=Structured datatype containing separate copies of all data of gs_hamk_in upon exit.
SOURCE
1284 subroutine copy_hamiltonian(gs_hamk_out,gs_hamk_in) 1285 1286 !Arguments ------------------------------------ 1287 type(gs_hamiltonian_type),intent(in),target :: gs_hamk_in 1288 type(gs_hamiltonian_type),intent(out),target :: gs_hamk_out 1289 1290 !Local variables------------------------------- 1291 integer :: tmp2i(5) 1292 #if defined HAVE_FC_ISO_C_BINDING 1293 type(C_PTR) :: ham_ptr 1294 #endif 1295 1296 ! ************************************************************************* 1297 1298 DBG_ENTER("COLL") 1299 1300 !@gs_hamiltonian_type 1301 1302 gs_hamk_out%dimekb1 = gs_hamk_in%dimekb1 1303 gs_hamk_out%dimekb2 = gs_hamk_in%dimekb2 1304 gs_hamk_out%dimekbq = gs_hamk_in%dimekbq 1305 gs_hamk_out%istwf_k = gs_hamk_in%istwf_k 1306 gs_hamk_out%istwf_kp = gs_hamk_in%istwf_kp 1307 gs_hamk_out%lmnmax = gs_hamk_in%lmnmax 1308 gs_hamk_out%matblk = gs_hamk_in%matblk 1309 gs_hamk_out%mgfft = gs_hamk_in%mgfft 1310 gs_hamk_out%mpsang = gs_hamk_in%mpsang 1311 gs_hamk_out%mpssoang = gs_hamk_in%mpssoang 1312 gs_hamk_out%natom = gs_hamk_in%natom 1313 gs_hamk_out%nfft = gs_hamk_in%nfft 1314 gs_hamk_out%npw_k = gs_hamk_in%npw_k 1315 gs_hamk_out%npw_kp = gs_hamk_in%npw_kp 1316 gs_hamk_out%npw_fft_k = gs_hamk_in%npw_fft_k 1317 gs_hamk_out%npw_fft_kp = gs_hamk_in%npw_fft_kp 1318 gs_hamk_out%nspinor = gs_hamk_in%nspinor 1319 gs_hamk_out%nsppol = gs_hamk_in%nsppol 1320 gs_hamk_out%ntypat = gs_hamk_in%ntypat 1321 gs_hamk_out%nvloc = gs_hamk_in%nvloc 1322 gs_hamk_out%n4 = gs_hamk_in%n4 1323 gs_hamk_out%n5 = gs_hamk_in%n5 1324 gs_hamk_out%n6 = gs_hamk_in%n6 1325 gs_hamk_out%gpu_option = gs_hamk_in%gpu_option 1326 gs_hamk_out%usecprj = gs_hamk_in%usecprj 1327 gs_hamk_out%usepaw = gs_hamk_in%usepaw 1328 gs_hamk_out%useylm = gs_hamk_in%useylm 1329 gs_hamk_out%ngfft = gs_hamk_in%ngfft 1330 gs_hamk_out%nloalg = gs_hamk_in%nloalg 1331 gs_hamk_out%ucvol = gs_hamk_in%ucvol 1332 gs_hamk_out%gmet = gs_hamk_in%gmet 1333 gs_hamk_out%gprimd = gs_hamk_in%gprimd 1334 gs_hamk_out%kpt_k = gs_hamk_in%kpt_k 1335 gs_hamk_out%kpt_kp = gs_hamk_in%kpt_kp 1336 1337 ABI_MALLOC(gs_hamk_out%atindx,(gs_hamk_out%natom)) 1338 gs_hamk_out%atindx = gs_hamk_in%atindx 1339 ABI_MALLOC(gs_hamk_out%atindx1,(gs_hamk_out%natom)) 1340 gs_hamk_out%atindx1 = gs_hamk_in%atindx1 1341 ABI_MALLOC(gs_hamk_out%dimcprj,(gs_hamk_out%natom*gs_hamk_out%usepaw)) 1342 if (gs_hamk_out%usepaw==1) gs_hamk_out%dimcprj = gs_hamk_in%dimcprj 1343 ABI_MALLOC(gs_hamk_out%typat,(gs_hamk_out%natom)) 1344 gs_hamk_out%typat = gs_hamk_in%typat 1345 ABI_MALLOC(gs_hamk_out%gbound_k,(2*gs_hamk_out%mgfft+8,2)) 1346 gs_hamk_out%gbound_k = gs_hamk_in%gbound_k 1347 ABI_MALLOC(gs_hamk_out%indlmn,(6,gs_hamk_out%lmnmax,gs_hamk_out%ntypat)) 1348 gs_hamk_out%indlmn = gs_hamk_in%indlmn 1349 ABI_MALLOC(gs_hamk_out%nattyp,(gs_hamk_out%ntypat)) 1350 gs_hamk_out%nattyp = gs_hamk_in%nattyp 1351 ABI_MALLOC(gs_hamk_out%nucdipmom,(3,gs_hamk_out%natom)) 1352 gs_hamk_out%nucdipmom = gs_hamk_in%nucdipmom 1353 ABI_MALLOC(gs_hamk_out%phkxred,(2,gs_hamk_out%natom)) 1354 gs_hamk_out%phkxred = gs_hamk_in%phkxred 1355 ABI_MALLOC(gs_hamk_out%ph1d,(2,3*(2*gs_hamk_out%mgfft+1)*gs_hamk_out%natom)) 1356 gs_hamk_out%ph1d = gs_hamk_in%ph1d 1357 ABI_MALLOC(gs_hamk_out%pspso,(gs_hamk_out%ntypat)) 1358 gs_hamk_out%pspso = gs_hamk_in%pspso 1359 tmp2i(1:5)=shape(gs_hamk_in%ekb_spin) 1360 ABI_MALLOC(gs_hamk_out%ekb_spin,(tmp2i(1),tmp2i(2),tmp2i(3),tmp2i(4),tmp2i(5))) 1361 gs_hamk_out%ekb_spin = gs_hamk_in%ekb_spin 1362 gs_hamk_out%ekb => gs_hamk_out%ekb_spin(:,:,:,:,1) 1363 tmp2i(1:2)=shape(gs_hamk_in%sij) 1364 ABI_MALLOC(gs_hamk_out%sij,(tmp2i(1),tmp2i(2))) 1365 gs_hamk_out%sij = gs_hamk_in%sij 1366 1367 if (associated(gs_hamk_in%gbound_kp,gs_hamk_in%gbound_k)) then 1368 gs_hamk_out%gbound_kp => gs_hamk_out%gbound_k 1369 else 1370 ABI_MALLOC(gs_hamk_out%gbound_kp,(2,gs_hamk_out%natom)) 1371 gs_hamk_out%gbound_kp = gs_hamk_in%gbound_kp 1372 end if 1373 if (associated(gs_hamk_in%phkpxred,gs_hamk_in%phkxred)) then 1374 gs_hamk_out%phkpxred => gs_hamk_out%phkxred 1375 else 1376 ABI_MALLOC(gs_hamk_out%phkpxred,(2,gs_hamk_out%natom)) 1377 gs_hamk_out%phkpxred = gs_hamk_in%phkpxred 1378 end if 1379 1380 call addr_copy(gs_hamk_in%xred,gs_hamk_out%xred) 1381 call addr_copy(gs_hamk_in%vectornd,gs_hamk_out%vectornd) 1382 call addr_copy(gs_hamk_in%vlocal,gs_hamk_out%vlocal) 1383 call addr_copy(gs_hamk_in%vxctaulocal,gs_hamk_out%vxctaulocal) 1384 call addr_copy(gs_hamk_in%kinpw_k,gs_hamk_out%kinpw_k) 1385 call addr_copy(gs_hamk_in%kinpw_kp,gs_hamk_out%kinpw_kp) 1386 call addr_copy(gs_hamk_in%kg_k,gs_hamk_out%kg_k) 1387 call addr_copy(gs_hamk_in%kg_kp,gs_hamk_out%kg_kp) 1388 call addr_copy(gs_hamk_in%kpg_k,gs_hamk_out%kpg_k) 1389 call addr_copy(gs_hamk_in%kpg_kp,gs_hamk_out%kpg_kp) 1390 call addr_copy(gs_hamk_in%ffnl_k,gs_hamk_out%ffnl_k) 1391 call addr_copy(gs_hamk_in%ffnl_kp,gs_hamk_out%ffnl_kp) 1392 call addr_copy(gs_hamk_in%ph3d_k,gs_hamk_out%ph3d_k) 1393 call addr_copy(gs_hamk_in%ph3d_kp,gs_hamk_out%ph3d_kp) 1394 1395 !For pointers to structured datatypes, have to copy the address 1396 !manually because there is no generic addr_copy function for that 1397 if (associated(gs_hamk_in%fockcommon)) then 1398 #if defined HAVE_FC_ISO_C_BINDING 1399 ham_ptr=c_loc(gs_hamk_in%fockcommon) 1400 call c_f_pointer(ham_ptr,gs_hamk_out%fockcommon) 1401 #else 1402 gs_hamk_out%fockcommon=transfer(gs_hamk_in%fockcommon,gs_hamk_out%fockcommon) 1403 #endif 1404 else 1405 nullify(gs_hamk_out%fockcommon) 1406 end if 1407 if (associated(gs_hamk_in%fockbz)) then 1408 #if defined HAVE_FC_ISO_C_BINDING 1409 ham_ptr=c_loc(gs_hamk_in%fockbz) 1410 call c_f_pointer(ham_ptr,gs_hamk_out%fockbz) 1411 #else 1412 gs_hamk_out%fockbz=transfer(gs_hamk_in%fockbz,gs_hamk_out%fockbz) 1413 #endif 1414 else 1415 nullify(gs_hamk_out%fockbz) 1416 end if 1417 if (associated(gs_hamk_in%fockACE_k)) then 1418 #if defined HAVE_FC_ISO_C_BINDING 1419 ham_ptr=c_loc(gs_hamk_in%fockACE_k) 1420 call c_f_pointer(ham_ptr,gs_hamk_out%fockACE_k) 1421 #else 1422 gs_hamk_out%fockACE_k=transfer(gs_hamk_in%fockACE_k,gs_hamk_out%fockACE_k) 1423 #endif 1424 else 1425 nullify(gs_hamk_out%fockACE_k) 1426 end if 1427 1428 DBG_EXIT("COLL") 1429 1430 end subroutine copy_hamiltonian
m_hamiltonian/destroy_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
destroy_hamiltonian
FUNCTION
Clean and destroy gs_hamiltonian_type datastructure
SIDE EFFECTS
Ham<gs_hamiltonian_type>=All dynamic memory defined in the structure is deallocated.
SOURCE
584 subroutine destroy_hamiltonian(Ham) 585 586 !Arguments ------------------------------------ 587 !scalars 588 class(gs_hamiltonian_type),intent(inout),target :: Ham 589 590 ! ************************************************************************* 591 592 DBG_ENTER("COLL") 593 594 !@gs_hamiltonian_type 595 596 ! Integer Pointers 597 if (associated(Ham%gbound_kp,Ham%gbound_k)) then 598 nullify(Ham%gbound_kp) 599 else if (associated(Ham%gbound_kp)) then 600 ABI_FREE(Ham%gbound_kp) 601 end if 602 603 ! Integer arrays 604 if(Ham%gpu_option == ABI_GPU_KOKKOS) then 605 #if defined HAVE_GPU && defined HAVE_YAKL 606 ABI_SFREE_MANAGED(Ham%atindx) 607 ABI_SFREE_MANAGED(Ham%atindx1) 608 ABI_SFREE_MANAGED(Ham%typat) 609 ABI_SFREE_MANAGED(Ham%indlmn) 610 ABI_SFREE_MANAGED(Ham%nattyp) 611 #endif 612 else 613 ABI_FREE(Ham%atindx) 614 ABI_FREE(Ham%atindx1) 615 ABI_FREE(Ham%typat) 616 ABI_FREE(Ham%indlmn) 617 ABI_FREE(Ham%nattyp) 618 end if 619 ABI_SFREE(Ham%gbound_k) 620 ABI_SFREE(Ham%pspso) 621 622 ABI_SFREE(Ham%dimcprj) 623 624 ! Real Pointers 625 if (associated(Ham%phkpxred,Ham%phkxred)) then 626 nullify(Ham%phkpxred) 627 else if (associated(Ham%phkpxred)) then 628 ABI_FREE(Ham%phkpxred) 629 end if 630 ABI_SFREE(Ham%phkxred) 631 if (associated(Ham%ekb)) nullify(Ham%ekb) 632 if (associated(Ham%vectornd)) nullify(Ham%vectornd) 633 if (associated(Ham%vlocal)) nullify(Ham%vlocal) 634 if (associated(Ham%vxctaulocal)) nullify(Ham%vxctaulocal) 635 if (associated(Ham%xred)) nullify(Ham%xred) 636 if (associated(Ham%kinpw_k)) nullify(Ham%kinpw_k) 637 if (associated(Ham%kinpw_kp)) nullify(Ham%kinpw_kp) 638 if (associated(Ham%kg_k)) nullify(Ham%kg_k) 639 if (associated(Ham%kg_kp)) nullify(Ham%kg_kp) 640 if (associated(Ham%kpg_k)) nullify(Ham%kpg_k) 641 if (associated(Ham%kpg_kp)) nullify(Ham%kpg_kp) 642 if (associated(Ham%ffnl_k)) nullify(Ham%ffnl_k) 643 if (associated(Ham%ffnl_kp)) nullify(Ham%ffnl_kp) 644 if (associated(Ham%ph3d_k)) nullify(Ham%ph3d_k) 645 if (associated(Ham%ph3d_kp)) nullify(Ham%ph3d_kp) 646 647 ! Real arrays 648 ABI_SFREE(Ham%ekb_spin) 649 ABI_SFREE(Ham%sij) 650 ABI_SFREE(Ham%nucdipmom) 651 if(Ham%gpu_option == ABI_GPU_KOKKOS) then 652 #if defined HAVE_GPU && defined HAVE_YAKL 653 ABI_SFREE_MANAGED(Ham%ph1d) 654 #endif 655 else 656 ABI_FREE(Ham%ph1d) 657 end if 658 659 ! Structured datatype pointers 660 if (associated(Ham%fockcommon)) nullify(Ham%fockcommon) 661 if (associated(Ham%fockACE_k)) nullify(Ham%fockACE_k) 662 if (associated(Ham%fockbz)) nullify(Ham%fockbz) 663 #if defined HAVE_GPU_CUDA 664 if(Ham%gpu_option==ABI_GPU_LEGACY .or. Ham%gpu_option==ABI_GPU_KOKKOS) then 665 call gpu_finalize_ham_data() 666 end if 667 #endif 668 669 DBG_EXIT("COLL") 670 671 end subroutine destroy_hamiltonian
m_hamiltonian/destroy_rf_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
destroy_rf_hamiltonian
FUNCTION
Clean and destroy rf_hamiltonian_type datastructure
SIDE EFFECTS
rf_Ham<rf_hamiltonian_type>=All dynamic memory defined in the structure is deallocated.
SOURCE
1526 subroutine destroy_rf_hamiltonian(rf_Ham) 1527 1528 !Arguments ------------------------------------ 1529 !scalars 1530 class(rf_hamiltonian_type),intent(inout) :: rf_Ham 1531 1532 ! ************************************************************************* 1533 1534 DBG_ENTER("COLL") 1535 1536 !@rf_hamiltonian_type 1537 1538 ! Real arrays 1539 ABI_SFREE(rf_Ham%e1kbfr_spin) 1540 ABI_SFREE(rf_Ham%e1kbsc_spin) 1541 1542 ! Real pointers 1543 if (associated(rf_Ham%dkinpw_k)) nullify(rf_Ham%dkinpw_k) 1544 if (associated(rf_Ham%dkinpw_kp)) nullify(rf_Ham%dkinpw_kp) 1545 if (associated(rf_Ham%ddkinpw_k)) nullify(rf_Ham%ddkinpw_k) 1546 if (associated(rf_Ham%ddkinpw_kp)) nullify(rf_Ham%ddkinpw_kp) 1547 if (associated(rf_Ham%vectornd)) nullify(rf_Ham%vectornd) 1548 if (associated(rf_Ham%vlocal1)) nullify(rf_Ham%vlocal1) 1549 if (associated(rf_Ham%vxctaulocal)) nullify(rf_Ham%vxctaulocal) 1550 if (associated(rf_Ham%e1kbfr)) nullify(rf_Ham%e1kbfr) 1551 if (associated(rf_Ham%e1kbsc)) nullify(rf_Ham%e1kbsc) 1552 1553 DBG_EXIT("COLL") 1554 1555 end subroutine destroy_rf_hamiltonian
m_hamiltonian/gs_hamiltonian_type [ Types ]
[ Top ] [ m_hamiltonian ] [ Types ]
NAME
gs_hamiltonian_type
FUNCTION
This datastructure contains the information about one Hamiltonian, needed in the "getghc" routine, that apply the Hamiltonian on a wavefunction. The Hamiltonian is expressed in reciprocal space: H_k^prime,k = exp(-i.k^prime.r^prime) H exp(i.k.r) In most cases k = k^prime and the k^prime objects are simply pointers to k objects.
SOURCE
96 type,public :: gs_hamiltonian_type 97 98 ! ===== Integer scalars 99 100 integer :: dimekb1 101 ! First dimension of Ekb 102 ! Same as psps%dimekb 103 ! ->Norm conserving : Max. number of Kleinman-Bylander energies 104 ! for each atom type 105 ! dimekb1=lnmax 106 ! ->PAW : Max. number of Dij coefficients connecting projectors 107 ! for each atom 108 ! dimekb1=cplex_dij*lmnmax*(lmnmax+1)/2 109 110 integer :: dimekb2 111 ! Second dimension of Ekb 112 ! ->Norm conserving psps: dimekb2=ntypat 113 ! ->PAW : dimekb2=natom 114 115 integer :: dimekbq 116 ! Fourth dimension of Ekb 117 ! 2 if Ekb factors contain a exp(-iqR) phase, 1 otherwise 118 119 integer :: istwf_k 120 ! option parameter that describes the storage of wfs at k 121 122 integer :: istwf_kp 123 ! option parameter that describes the storage of wfs at k^prime 124 125 integer :: lmnmax 126 ! Maximum number of different l,m,n components over all types of psps. 127 ! same as dtset%lmnmax 128 129 integer :: matblk 130 ! dimension of the array ph3d 131 132 integer :: mgfft 133 ! maximum size for 1D FFTs (same as dtset%mgfft) 134 135 integer :: mpsang 136 ! Highest angular momentum of non-local projectors over all type of psps. 137 ! shifted by 1 : for all local psps, mpsang=0; for largest s, mpsang=1, 138 ! for largest p, mpsang=2; for largest d, mpsang=3; for largest f, mpsang=4 139 ! This gives also the number of non-local "channels" 140 ! same as psps%mpsang 141 142 integer :: mpssoang 143 ! Maximum number of channels, including those for treating the spin-orbit coupling 144 ! For NC pseudopotentials only: 145 ! when mpspso=1, mpssoang=mpsang 146 ! when mpspso=2, mpssoang=2*mpsang-1 147 ! For PAW: same as mpsang 148 ! same as psps%mpssoang 149 150 integer :: natom 151 ! The number of atoms for this dataset; same as dtset%natom 152 153 integer :: nfft 154 ! number of FFT grid points same as dtset%nfft 155 156 integer :: npw_k 157 ! number of plane waves at k 158 ! In case of band-FFT parallelism, npw_k is the number of plane waves 159 ! processed by current proc 160 161 integer :: npw_fft_k 162 ! number of plane waves at k used to apply Hamiltonian when band-FFT 163 ! parallelism is activated (i.e. data are distributed in the "FFT" configuration) 164 165 integer :: npw_kp 166 ! number of plane waves at k^prime 167 ! In case of band-FFT parallelism, npw_kp is the number of plane waves 168 ! processed by current proc 169 170 integer :: npw_fft_kp 171 ! number of plane waves at k^prime used to apply Hamiltonian when band-FFT 172 ! parallelism is activated (i.e. data are distributed in the "FFT" configuration) 173 174 integer :: nspinor 175 ! Number of spinorial components 176 177 integer :: nsppol 178 ! Total number of spin components (1=non-polarized, 2=polarized) 179 180 integer :: ntypat 181 ! Number of types of pseudopotentials same as dtset%ntypat 182 183 integer :: nvloc 184 ! Number of components of vloc 185 ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4 186 187 integer :: n4,n5,n6 188 ! same as ngfft(4:6) 189 190 integer :: gpu_option 191 ! Governs the choice of the GPU implementation: 192 ! = 0 ==> do not use GPU 193 ! > 0 ==> see defs_basis.F90 to have the list of possible GPU implementations 194 195 integer :: usecprj 196 ! usecprj= 1 if cprj projected WF are stored in memory 197 ! = 0 if they are to be computed on the fly 198 199 integer :: usepaw 200 ! if usepaw=0 , use norm-conserving psps part of the code 201 ! is usepaw=1 , use paw part of the code 202 203 integer :: useylm 204 ! governs the way the nonlocal operator is to be applied: 205 ! 1=using Ylm, 0=using Legendre polynomials 206 207 ! ===== Integer arrays 208 209 #if defined HAVE_GPU && defined HAVE_YAKL 210 integer(c_int32_t), ABI_CONTIGUOUS pointer :: atindx(:) => null() 211 #else 212 integer, allocatable :: atindx(:) 213 #endif 214 ! atindx(natom) 215 ! index table for atoms (see gstate.f) 216 217 #if defined HAVE_GPU && defined HAVE_YAKL 218 integer(c_int32_t), ABI_CONTIGUOUS pointer :: atindx1(:) => null() 219 #else 220 integer, allocatable :: atindx1(:) 221 #endif 222 ! atindx1(natom) 223 ! index table for atoms, inverse of atindx (see gstate.f) 224 225 integer, allocatable :: dimcprj(:) 226 ! dimcprj(natom*usepaw)=dimensions of array cprj 227 ! dimcprj(ia)=cprj(ia,:)%nlmn 228 ! atoms are ordered by atom-type 229 230 integer, allocatable :: gbound_k(:,:) 231 ! gbound_k(2*mgfft+8,2) 232 ! G sphere boundary, for each plane wave at k 233 234 #if defined HAVE_GPU && defined HAVE_YAKL 235 integer(c_int32_t), ABI_CONTIGUOUS pointer :: indlmn(:,:,:) => null() 236 #else 237 integer(c_int32_t), allocatable :: indlmn(:,:,:) 238 #endif 239 ! indlmn(6,lmnmax,ntypat) 240 ! For each type of psp, 241 ! array giving l,m,n,lm,ln,spin for i=ln (if useylm=0) 242 ! or i=lmn (if useylm=1) 243 244 #if defined HAVE_GPU && defined HAVE_YAKL 245 integer(c_int32_t), ABI_CONTIGUOUS pointer :: nattyp(:) => null() 246 #else 247 integer, allocatable :: nattyp(:) 248 #endif 249 ! nattyp(ntypat) 250 ! # of atoms of each type 251 252 integer :: ngfft(18) 253 ! ngfft(1:3)=integer fft box dimensions 254 ! ngfft(4:6)=integer fft box dimensions, might be augmented for CPU speed 255 ! ngfft(7)=fftalg 256 ! ngfft(8)=fftalg 257 258 integer :: nloalg(3) 259 ! governs the choice of the algorithm for non-local operator same as dtset%nloalg 260 261 integer, allocatable :: pspso(:) 262 ! pspso(ntypat) 263 ! For each type of psp, 1 if no spin-orbit component is taken 264 ! into account, 2 if a spin-orbit component is used 265 ! Revelant for NC-psps and PAW 266 267 #if defined HAVE_GPU && defined HAVE_YAKL 268 integer(c_int32_t), ABI_CONTIGUOUS pointer :: typat(:) => null() 269 #else 270 integer, allocatable :: typat(:) 271 #endif 272 ! typat(natom) 273 ! type of each atom 274 275 ! integer, allocatable :: indpw_k(:,:) 276 ! indpw_k(4,npw_fft_k) 277 ! array which gives fft box index for given basis sphere 278 279 ! Integer pointers 280 281 integer, ABI_CONTIGUOUS pointer :: gbound_kp(:,:) => null() 282 ! gbound_kp(2*mgfft+8,2) 283 ! G sphere boundary, for each plane wave at k^prime 284 285 #if defined HAVE_GPU && defined HAVE_YAKL 286 integer(int32), ABI_CONTIGUOUS pointer :: kg_k(:,:) => null() 287 #else 288 integer, pointer :: kg_k(:,:) => null() 289 #endif 290 ! kg_k(3,npw_fft_k) 291 ! G vector coordinates with respect to reciprocal lattice translations 292 ! at k 293 294 integer, pointer :: kg_kp(:,:) => null() 295 ! kg_kp(3,npw_fft_kp) 296 ! G vector coordinates with respect to reciprocal lattice translations 297 ! at k^prime 298 299 ! ===== Real scalars 300 301 real(dp) :: ucvol 302 ! unit cell volume (Bohr**3) 303 304 ! ===== Real arrays 305 306 real(dp), allocatable :: ekb_spin(:,:,:,:,:) 307 ! ekb_spin(dimekb1,dimekb2,nspinor**2,dimekbq,my_nsppol) 308 ! Contains the values of ekb array for all spins treated by current process 309 ! See ekb description ; ekb is pointer to ekb_spin(:,:,:,:,my_isppol) 310 311 real(dp), allocatable :: sij(:,:) 312 ! sij(dimekb1,ntypat*usepaw) = overlap matrix for paw calculation 313 314 real(dp) :: gmet(3,3) 315 ! reciprocal space metric tensor in Bohr**-2 316 317 real(dp) :: gprimd(3,3) 318 ! dimensional reciprocal space primitive translations (Bohr^-1) 319 320 real(dp) :: kpt_k(3) 321 ! dimensionless k point coordinates wrt reciprocal lattice vectors 322 323 real(dp) :: kpt_kp(3) 324 ! dimensionless k^prime point coordinates wrt reciprocal lattice vectors 325 326 real(dp), allocatable :: nucdipmom(:,:) 327 ! nucdipmom(3,natom) 328 ! nuclear dipole moments at each atomic position 329 330 #if defined HAVE_GPU && defined HAVE_YAKL 331 real(c_double), ABI_CONTIGUOUS pointer :: ph1d(:,:) => null() 332 #else 333 real(dp), allocatable :: ph1d(:,:) 334 #endif 335 ! ph1d(2,3*(2*mgfft+1)*natom) 336 ! 1-dim phase arrays for structure factor (see getph.f). 337 338 real(dp), allocatable :: phkxred(:,:) 339 ! phkxred(2,natom) 340 ! phase factors exp(2 pi k.xred) at k 341 342 ! ===== Real pointers 343 344 real(dp), ABI_CONTIGUOUS pointer :: ekb(:,:,:,:) => null() 345 ! ekb(dimekb1,dimekb2,nspinor**2,dimekbq) 346 ! ->Norm conserving : (Real) Kleinman-Bylander energies (hartree) 347 ! for number of basis functions (l,n) (lnmax) 348 ! and number of atom types (ntypat) 349 ! dimekb1=lnmax ; dimekb2=ntypat ; dimekbq=1 350 ! ->PAW : (Real, symmetric) Frozen part of Dij coefficients 351 ! to connect projectors 352 ! for number of basis functions (l,m,n) (lmnmax) 353 ! and number of atom (natom) 354 ! dimekb1=lmnmax*(lmnmax+1)/2 ; dimekb2=natom ; dimekbq=1 355 ! ekb is spin dependent in the case of PAW calculations. 356 ! For each spin component, ekb points to ekb_spin(:,:,:,:,my_isppol) 357 ! dimekbq=2 if Ekb factors contain a exp(-iqR) phase, dimekbq=1 otherwise 358 ! About the non-local factors symmetry: 359 ! - The lower triangular part of the Dij matrix can be deduced from the upper one 360 ! with the following relation: D^s2s1_ji = (D^s1s2_ij)^* 361 ! where s1,s2 are spinor components 362 363 real(dp), pointer :: ffnl_k(:,:,:,:) => null() 364 ! ffnl_k(npw_fft_k,2,dimffnl_k,ntypat) 365 ! nonlocal form factors at k 366 367 real(dp), pointer :: ffnl_kp(:,:,:,:) => null() 368 ! ffnl_kp(npw_fft_kp,2,dimffnl_kp,ntypat) 369 ! nonlocal form factors at k_prime 370 371 real(dp), pointer :: kinpw_k(:) => null() 372 ! kinpw_k(npw_fft_k) 373 ! (modified) kinetic energy for each plane wave at k 374 ! CAVEAT: In band mode, this array is NOT EQUIVALENT to kinpw(npw_k) 375 376 real(dp), pointer :: kinpw_kp(:) => null() 377 ! kinpw_kp(npw_fft_kp) 378 ! (modified) kinetic energy for each plane wave at k^prime 379 380 real(dp), pointer :: kpg_k(:,:) => null() 381 ! kpg_k(3,npw_fft_k) 382 ! k+G vector coordinates at k 383 384 real(dp), pointer :: kpg_kp(:,:) => null() 385 ! kpg_kp(3,npw_fft_kp) 386 ! k^prime+G vector coordinates at k^prime 387 388 real(dp), ABI_CONTIGUOUS pointer :: phkpxred(:,:) => null() 389 ! phkpxred(2,natom) 390 ! phase factors exp(2 pi k^prime.xred) at k^prime 391 392 real(dp), pointer :: ph3d_k(:,:,:) => null() 393 ! ph3d_k(2,npw_fft_k,matblk) 394 ! 3-dim structure factors, for each atom and plane wave at k 395 396 real(dp), pointer :: ph3d_kp(:,:,:) => null() 397 ! ph3d_kp(2,npw_fft_kp,matblk) 398 ! 3-dim structure factors, for each atom and plane wave at k^prime 399 400 real(dp), pointer :: vectornd(:,:,:,:,:) => null() 401 ! vectornd(n4,n5,n6,nvloc,3) 402 ! vector potential of nuclear magnetic dipoles 403 ! in real space, on the augmented fft grid 404 405 real(dp), pointer :: vlocal(:,:,:,:) => null() 406 ! vlocal(n4,n5,n6,nvloc) 407 ! local potential in real space, on the augmented fft grid 408 409 real(dp), pointer :: vxctaulocal(:,:,:,:,:) => null() 410 ! vxctaulocal(n4,n5,n6,nvloc,4) 411 ! derivative of XC energy density with respect to kinetic energy density, 412 ! in real space, on the augmented fft grid 413 414 real(dp), ABI_CONTIGUOUS pointer :: xred(:,:) => null() 415 ! xred(3,natom) 416 ! reduced coordinates of atoms (dimensionless) 417 418 ! ===== Structured datatype pointers 419 420 type(fock_common_type), pointer :: fockcommon => null() 421 ! common quantities needed to calculate Fock exact exchange 422 423 type(fock_BZ_type), pointer :: fockbz => null() 424 ! total brillouin zone quantities needed to calculate Fock exact exchange 425 426 type(fock_ACE_type), pointer :: fockACE_k => null() 427 ! ACE quantities needed to calculate Fock exact exchange in the ACE context 428 429 contains 430 procedure :: free => destroy_hamiltonian 431 ! Free the memory in the GS Hamiltonian 432 433 procedure :: load_spin => load_spin_hamiltonian 434 ! Setup of the spin-dependent part of the GS Hamiltonian 435 436 procedure :: load_k => load_k_hamiltonian 437 ! Setup of the k-dependent part of the GS Hamiltonian 438 439 procedure :: load_kprime => load_kprime_hamiltonian 440 ! Setup of the k^prime-dependent part of the GS Hamiltonian 441 442 !procedure :: copy => copy_hamiltonian 443 444 end type gs_hamiltonian_type 445 446 public :: init_hamiltonian ! Initialize the GS Hamiltonian 447 public :: copy_hamiltonian ! Copy the object
m_hamiltonian/init_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
init_hamiltonian
FUNCTION
Creation method for the gs_hamiltonian_type structure. It allocates memory and initializes all quantities that do not depend on the k-point or spin.
INPUTS
[comm_atom]=optional, MPI communicator over atoms [fock <type(fock_type)>]= common quantities to calculate Fock exact exchange natom=Number of atoms in the unit cell. nfft=Number of FFT grid points (for this processors). nspinor=Number of spinorial components nsppol=1 for unpolarized, 2 for spin-polarized nspden=Number of spin density components. mgfft=Maximum size for 1D FFTs i.e., MAXVAL(ngfft(1:3)) [mpi_atmtab(:)]=optional, indexes of the atoms treated by current proc [mpi_spintab(2)]=optional, flags defining the spin(s) treated be current process: mpi_spintab(1)=1 if non-polarized or spin-up treated mpi_spintab(2)=1 if polarized and spin-dn treated psps<pseudopotential_type>=structure datatype gathering data on the pseudopotentials. [electronpositron<electronpositron_type>]=Structured datatype storing data for the electron-positron two-component DFT (optional). ngfft(18)=integer array with FFT box dimensions and other information on FFTs, for the FINE rectangular grid. nloalg(3)=governs the choice of the algorithm for non-local operator [nucdipmom(3,natom)]= (optional) array of nuclear dipole moments at atomic sites [ph1d(2,3*(2*mgfft+1)*natom)]=1-dimensions phase arrays for structure factor (see getph.f). Optional, recalculated inside the routine if not present in input. rprimd(3,3)=Direct lattice vectors in Bohr. typat(natom)=Type of each atom. [usecprj]=flag use only for PAW; 1 if cprj datastructure is allocated [gpu_option] = GPU implementation to use, i.e. cuda, openMP, ... (0=not using GPU) xred(3,natom)=Reduced coordinates of the atoms. pawtab(ntypat*psps%usepaw)<pawtab_type>=PAW TABulated data initialized at start. [paw_ij(:) <type(paw_ij_type)>]=optional, paw arrays given on (i,j) channels
SIDE EFFECTS
Ham<gs_hamiltonian_type>=Structured datatype almost completely initialized: * Basic variables and dimensions are transfered to the structure. * All pointers are allocated with correct dimensions. * Quantities that do not depend on the k-point or spin are initialized.
SOURCE
721 subroutine init_hamiltonian(ham,Psps,pawtab,nspinor,nsppol,nspden,natom,typat,& 722 & xred,nfft,mgfft,ngfft,rprimd,nloalg,& 723 & ph1d,usecprj,comm_atom,mpi_atmtab,mpi_spintab,paw_ij,& ! optional 724 & electronpositron,fock,nucdipmom,gpu_option) ! optional 725 726 !Arguments ------------------------------------ 727 !scalars 728 integer,intent(in) :: nfft,natom,nspinor,nsppol,nspden,mgfft 729 integer,optional,intent(in) :: comm_atom,usecprj,gpu_option 730 class(gs_hamiltonian_type),intent(inout),target :: ham 731 type(electronpositron_type),optional,pointer :: electronpositron 732 type(fock_type),optional,pointer :: fock 733 type(pseudopotential_type),intent(in) :: psps 734 !arrays 735 integer,intent(in) :: ngfft(18),nloalg(3),typat(natom) 736 integer,optional,intent(in) :: mpi_atmtab(:),mpi_spintab(2) 737 real(dp),intent(in) :: rprimd(3,3) 738 real(dp),intent(in),target :: xred(3,natom) 739 real(dp),optional,intent(in) :: nucdipmom(3,natom),ph1d(2,3*(2*mgfft+1)*natom) 740 type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw) 741 type(paw_ij_type),optional,intent(in) :: paw_ij(:) 742 743 !Local variables------------------------------- 744 !scalars 745 integer :: my_comm_atom,my_nsppol,itypat,iat,ilmn,indx,isp,cplex_dij,jsp,l_gpu_option 746 real(dp) :: ucvol 747 !arrays 748 integer :: my_spintab(2) 749 real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3) 750 real(dp),allocatable,target :: ekb_tmp(:,:,:,:) 751 752 ! ************************************************************************* 753 754 DBG_ENTER("COLL") 755 756 !@gs_hamiltonian_type 757 758 !Manage optional parameters 759 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 760 my_spintab=0;my_spintab(1:nsppol)=1;if (present(mpi_spintab)) my_spintab(1:2)=mpi_spintab(1:2) 761 my_nsppol=count(my_spintab==1) 762 l_gpu_option=ABI_GPU_DISABLED; if(present(gpu_option)) l_gpu_option=gpu_option 763 764 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 765 766 ABI_CHECK(mgfft==MAXVAL(ngfft(1:3)),"Wrong mgfft") 767 768 !Allocate the arrays of the Hamiltonian whose dimensions do not depend on k 769 if(l_gpu_option == ABI_GPU_KOKKOS) then 770 #if defined HAVE_GPU && defined HAVE_YAKL 771 ABI_MALLOC_MANAGED(ham%atindx,(/natom/)) 772 ABI_MALLOC_MANAGED(ham%atindx1,(/natom/)) 773 ABI_MALLOC_MANAGED(ham%typat,(/natom/)) 774 ABI_MALLOC_MANAGED(ham%indlmn,(/6,psps%lmnmax,psps%ntypat/)) 775 ABI_MALLOC_MANAGED(ham%nattyp,(/psps%ntypat/)) 776 ABI_MALLOC_MANAGED(ham%ph1d,(/2,3*(2*mgfft+1)*natom/)) 777 #endif 778 else 779 ABI_MALLOC(ham%atindx,(natom)) 780 ABI_MALLOC(ham%atindx1,(natom)) 781 ABI_MALLOC(ham%typat,(natom)) 782 ABI_MALLOC(ham%indlmn,(6,psps%lmnmax,psps%ntypat)) 783 ABI_MALLOC(ham%nattyp,(psps%ntypat)) 784 ABI_MALLOC(ham%ph1d,(2,3*(2*mgfft+1)*natom)) 785 end if 786 787 ABI_MALLOC(ham%pspso,(psps%ntypat)) 788 ABI_MALLOC(ham%nucdipmom,(3,natom)) 789 790 !Initialize most of the Hamiltonian 791 indx=1 792 do itypat=1,psps%ntypat 793 ham%nattyp(itypat)=0 794 do iat=1,natom 795 if (typat(iat)==itypat) then 796 ham%atindx (iat )=indx 797 ham%atindx1(indx)=iat 798 indx=indx+1 799 ham%nattyp(itypat)=ham%nattyp(itypat)+1 800 end if 801 end do 802 end do 803 804 ham%gmet(:,:) =gmet(:,:) 805 ham%gprimd(:,:)=gprimd(:,:) 806 ham%indlmn(:,:,:)=psps%indlmn(:,:,:) 807 ham%lmnmax =psps%lmnmax 808 ham%mgfft =mgfft 809 ham%mpsang =psps%mpsang 810 ham%mpssoang =psps%mpssoang 811 ham%natom =natom 812 ham%nfft =nfft 813 ham%ngfft(:) =ngfft(:) 814 ham%nloalg(:) =nloalg(:) 815 ham%matblk=min(NLO_MINCAT,maxval(ham%nattyp)); if (nloalg(2)>0) ham%matblk=natom 816 ham%nsppol =nsppol 817 ham%nspinor =nspinor 818 ham%ntypat =psps%ntypat 819 ham%typat =typat(1:natom) 820 ham%nvloc=1; if(nspden==4)ham%nvloc=4 821 ham%n4 =ngfft(4) 822 ham%n5 =ngfft(5) 823 ham%n6 =ngfft(6) 824 ham%usepaw =psps%usepaw 825 ham%ucvol =ucvol 826 ham%useylm =psps%useylm 827 ham%gpu_option=ABI_GPU_DISABLED ; if(PRESENT(gpu_option)) ham%gpu_option=gpu_option 828 829 ham%pspso(:) =psps%pspso(1:psps%ntypat) 830 if (psps%usepaw==1) then 831 do itypat=1,psps%ntypat 832 ham%pspso(itypat)=1+pawtab(itypat)%usespnorb 833 end do 834 end if 835 836 if (present(nucdipmom)) then 837 ham%nucdipmom(:,:) = nucdipmom(:,:) 838 else 839 ham%nucdipmom(:,:) = zero 840 end if 841 842 ham%xred => xred 843 844 if (present(fock)) then 845 if (associated(fock)) then 846 ham%fockcommon => fock%fock_common 847 if (fock%fock_common%use_ACE==0) ham%fockbz=>fock%fock_BZ 848 end if 849 end if 850 851 if (present(ph1d)) then 852 ham%ph1d(:,:) = ph1d(:,:) 853 else ! Recalculate structure factor phases 854 call getph(ham%atindx,natom,ngfft(1),ngfft(2),ngfft(3),ham%ph1d,xred) 855 end if 856 857 if (ham%usepaw==1) then 858 ham%usecprj=0;if (present(usecprj)) ham%usecprj=usecprj 859 ABI_MALLOC(ham%dimcprj,(natom)) 860 !Be carefull cprj are ordered by atom type (used in non-local operator) 861 call pawcprj_getdim(ham%dimcprj,natom,ham%nattyp,ham%ntypat,ham%typat,pawtab,'O') 862 else 863 ham%usecprj=0 864 ABI_MALLOC(ham%dimcprj,(0)) 865 end if 866 867 ! =========================== 868 ! ==== Non-local factors ==== 869 ! =========================== 870 871 872 if (ham%usepaw==0) then ! Norm-conserving: use constant Kleimann-Bylander energies. 873 ham%dimekb1=psps%dimekb 874 ham%dimekb2=psps%ntypat 875 ham%dimekbq=1 876 ABI_MALLOC(ham%ekb_spin,(psps%dimekb,psps%ntypat,nspinor**2,1,1)) 877 ham%ekb => ham%ekb_spin(:,:,:,:,1) 878 ABI_MALLOC(ham%sij,(0,0)) 879 ham%ekb(:,:,1,1)=psps%ekb(:,:) 880 if (nspinor==2) then 881 ham%ekb(:,:,2,1)=psps%ekb(:,:) 882 ham%ekb(:,:,3:4,1)=zero 883 end if 884 if (PRESENT(electronpositron)) then 885 if (electronpositron_calctype(electronpositron)==1) ham%ekb(:,:,:,:)=-ham%ekb(:,:,:,:) 886 end if 887 888 ! Update enl on GPU (will do it later for PAW) 889 #if defined HAVE_GPU_CUDA 890 if (ham%gpu_option==ABI_GPU_LEGACY .or. ham%gpu_option==ABI_GPU_KOKKOS) then 891 call gpu_update_ham_data(& 892 & ham%ekb(:,:,:,1), INT(size(ham%ekb), c_int64_t), & 893 & ham%sij, INT(size(ham%sij), c_int64_t), & 894 & ham%gprimd, INT(size(ham%gprimd),c_int64_t)) 895 end if 896 #endif 897 898 else ! PAW: store overlap coefficients (spin non dependent) and Dij coefficients (spin dependent) 899 cplex_dij=1 900 if (present(paw_ij)) then 901 if (size(paw_ij)>0) cplex_dij=paw_ij(1)%cplex_dij 902 end if 903 if ((nspinor==2).or.any(abs(ham%nucdipmom)>tol8)) cplex_dij=2 904 ham%dimekb1=psps%dimekb*cplex_dij 905 ham%dimekb2=natom 906 ham%dimekbq=1 907 if (present(paw_ij)) then 908 if (size(paw_ij)>0) ham%dimekbq=paw_ij(1)%qphase 909 end if 910 ABI_MALLOC(ham%sij,(ham%dimekb1,psps%ntypat)) 911 do itypat=1,psps%ntypat 912 if (cplex_dij==1) then 913 ham%sij(1:pawtab(itypat)%lmn2_size,itypat)=pawtab(itypat)%sij(:) 914 else 915 do ilmn=1,pawtab(itypat)%lmn2_size 916 ham%sij(2*ilmn-1,itypat)=pawtab(itypat)%sij(ilmn) 917 ham%sij(2*ilmn ,itypat)=zero 918 end do 919 end if 920 if (cplex_dij*pawtab(itypat)%lmn2_size<ham%dimekb1) then 921 ham%sij(cplex_dij*pawtab(itypat)%lmn2_size+1:ham%dimekb1,itypat)=zero 922 end if 923 end do 924 !We preload here PAW non-local factors in order to avoid a communication over atoms 925 ! inside the loop over spins. 926 ABI_MALLOC(ham%ekb_spin,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq,my_nsppol)) 927 ham%ekb_spin=zero 928 if (present(paw_ij)) then 929 if (my_nsppol<ham%nsppol) then 930 ABI_MALLOC(ekb_tmp,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq)) 931 end if 932 jsp=0 933 do isp=1,ham%nsppol 934 if (my_spintab(isp)==1) then 935 jsp=jsp+1 ; ham%ekb => ham%ekb_spin(:,:,:,:,jsp) 936 else 937 ham%ekb => ekb_tmp 938 end if 939 if (present(mpi_atmtab)) then 940 call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom,mpi_atmtab=mpi_atmtab) 941 else 942 call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom) 943 end if 944 end do 945 if (my_nsppol<ham%nsppol) then 946 ABI_FREE(ekb_tmp) 947 end if 948 end if 949 nullify(ham%ekb) 950 end if 951 952 DBG_EXIT("COLL") 953 954 end subroutine init_hamiltonian
m_hamiltonian/init_rf_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
init_rf_hamiltonian
FUNCTION
Creation method for the rf_hamiltonian_type structure. It allocates memory and initializes all quantities that do not depend on the k-point or spin.
INPUTS
[comm_atom]=optional, MPI communicator over atoms cplex_paw=1 if all on-site PAW quantities are real (GS), 2 if they are complex (RF) gs_Ham<gs_hamiltonian_type>=Structured datatype containing data for ground-state Hamiltonian at (k+q) [has_e1kbsc]=optional, true if rf_Ham%e1kbsc has to be initialized. e1kbsc contains the self-consistent 1st-order PAW Dij coefficients (depending on VHxc^(1)) ipert=index of perturbation [mpi_atmtab(:)]=optional, indexes of the atoms treated by current proc [mpi_spintab(2)]=optional, flags defining the spin(s) treated be current process: mpi_spintab(1)=1 if non-polarized or spin-up treated mpi_spintab(2)=1 if polarized and spin-dn treated [paw_ij1(:)<paw_ij_type>]=Various 1st-order arrays given on (i,j) (partial waves) channels (paw_ij1%dij and paw_ij1%difr only used here).
SIDE EFFECTS
rf_Ham<rf_hamiltonian_type>=Structured datatype almost completely initialized: * Basic variables and dimensions are transfered to the structure. * All pointers are allocated with correct dimensions. * Quantities that do not depend on the k-point or spin are initialized.
SOURCE
1590 subroutine init_rf_hamiltonian(cplex,gs_Ham,ipert,rf_Ham,& 1591 & comm_atom,mpi_atmtab,mpi_spintab,paw_ij1,has_e1kbsc) ! optional arguments 1592 1593 !Arguments ------------------------------------ 1594 !scalars 1595 integer,intent(in) :: cplex,ipert 1596 integer,intent(in),optional :: comm_atom 1597 logical,intent(in),optional :: has_e1kbsc 1598 type(gs_hamiltonian_type),intent(in) :: gs_Ham 1599 type(rf_hamiltonian_type),intent(inout),target :: rf_Ham 1600 !arrays 1601 integer,optional,intent(in) :: mpi_atmtab(:),mpi_spintab(2) 1602 type(paw_ij_type),optional,intent(in) :: paw_ij1(:) 1603 1604 !Local variables------------------------------- 1605 !scalars 1606 integer :: cplex_dij1,isp,jsp,my_comm_atom,my_nsppol 1607 logical :: has_e1kbsc_ 1608 !arrays 1609 integer :: my_spintab(2) 1610 real(dp),allocatable,target :: e1kb_tmp(:,:,:,:) 1611 1612 ! ************************************************************************* 1613 1614 DBG_ENTER("COLL") 1615 1616 !@rf_hamiltonian_type 1617 1618 !Manage optional parameters 1619 has_e1kbsc_=.false.;if (present(has_e1kbsc)) has_e1kbsc_=has_e1kbsc 1620 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 1621 my_spintab=0;my_spintab(1:gs_Ham%nsppol)=1;if(present(mpi_spintab)) my_spintab=mpi_spintab 1622 my_nsppol=count(my_spintab==1) 1623 1624 rf_Ham%cplex =cplex 1625 rf_Ham%n4 =gs_Ham%n4 1626 rf_Ham%n5 =gs_Ham%n5 1627 rf_Ham%n6 =gs_Ham%n6 1628 rf_Ham%nvloc =gs_Ham%nvloc 1629 rf_Ham%nsppol =gs_Ham%nsppol 1630 rf_Ham%nspinor =gs_Ham%nspinor 1631 1632 rf_Ham%dime1kb1=0 1633 rf_Ham%dime1kb2=gs_Ham%dimekb2 1634 if (gs_Ham%usepaw==1.and.ipert/=gs_Ham%natom+1.and.ipert/=gs_Ham%natom+10) then 1635 cplex_dij1=1;if ((gs_Ham%nspinor==2).or.any(abs(gs_Ham%nucdipmom)>tol8)) cplex_dij1=2 1636 rf_Ham%dime1kb1=cplex_dij1*(gs_Ham%lmnmax*(gs_Ham%lmnmax+1))/2 1637 end if 1638 1639 ! Allocate the arrays of the 1st-order Hamiltonian 1640 ! We preload here 1st-order non-local factors in order to avoid 1641 ! a communication over atoms inside the loop over spins. 1642 if (gs_Ham%usepaw==1.and.rf_Ham%dime1kb1>0) then 1643 if ((ipert>=1.and.ipert<=gs_Ham%natom).or.ipert==gs_Ham%natom+2.or.& 1644 ipert==gs_Ham%natom+3.or.ipert==gs_Ham%natom+4.or.ipert==gs_Ham%natom+11) then 1645 1646 ABI_MALLOC(rf_Ham%e1kbfr_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol)) 1647 rf_Ham%e1kbfr_spin=zero 1648 if (has_e1kbsc_) then 1649 ABI_MALLOC(rf_Ham%e1kbsc_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol)) 1650 rf_Ham%e1kbsc_spin=zero 1651 end if 1652 1653 if (present(paw_ij1)) then 1654 1655 if (my_nsppol<rf_Ham%nsppol) then 1656 ABI_MALLOC(e1kb_tmp,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex)) 1657 end if 1658 1659 ! === Frozen term 1660 jsp=0 1661 do isp=1,rf_Ham%nsppol 1662 if (my_spintab(isp)==1) then 1663 jsp=jsp+1 ; rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsp) 1664 else 1665 rf_Ham%e1kbfr => e1kb_tmp 1666 end if 1667 if (present(mpi_atmtab)) then 1668 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr,mpi_atmtab=mpi_atmtab) 1669 else 1670 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr) 1671 end if 1672 end do 1673 1674 ! === Self-consistent term 1675 if (has_e1kbsc_) then 1676 jsp=0 1677 do isp=1,rf_Ham%nsppol 1678 if (my_spintab(isp)==1) then 1679 jsp=jsp+1 ; rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsp) 1680 else 1681 rf_Ham%e1kbsc => e1kb_tmp 1682 end if 1683 if (present(mpi_atmtab)) then 1684 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc,mpi_atmtab=mpi_atmtab) 1685 else 1686 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc) 1687 end if 1688 end do 1689 end if 1690 1691 if (my_nsppol<rf_Ham%nsppol) then 1692 ABI_FREE(e1kb_tmp) 1693 end if 1694 1695 end if 1696 end if 1697 end if 1698 1699 if (.not.allocated(rf_Ham%e1kbfr_spin)) then 1700 ABI_MALLOC(rf_Ham%e1kbfr_spin,(0,0,0,0,0)) 1701 end if 1702 if (.not.allocated(rf_Ham%e1kbsc_spin)) then 1703 ABI_MALLOC(rf_Ham%e1kbsc_spin,(0,0,0,0,0)) 1704 end if 1705 nullify(rf_Ham%e1kbfr) 1706 nullify(rf_Ham%e1kbsc) 1707 1708 DBG_EXIT("COLL") 1709 1710 end subroutine init_rf_hamiltonian
m_hamiltonian/load_k_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
load_k_hamiltonian
FUNCTION
Setup of the k-dependent part of the Hamiltonian H_k_k^prime
INPUTS
[compute_gbound]=flag. if true, G sphere boundary is computed here [compute_ph3d]=flag. if true, 3D structure factors are computed here (only if nloalg(1)>0) [gbound_k]=G sphere boundary (not compatible with compute_gbound=TRUE) [ffnl_k]=nonlocal form factors on basis sphere [istwf_k]=parameter that describes the storage of wfs [kinpw_k]=(modified) kinetic energy for each plane wave [kg_k]=planewave reduced coordinates in basis sphere (g vectors) [kpg_k]=(k+g) vectors in reciprocal space [kpt_k]=k point coordinates [npw_k]=number of plane waves (processed by current proc when band-FFT parallelism is on) [npw_fft_k]=number of plane waves used to apply Hamiltonian (in the "FFT" configuration) [ph3d_k]=3-dim structure factors, for each atom and plane wave
SIDE EFFECTS
ham<gs_hamiltonian_type>=structured datatype completed with k-dependent quantitites. Quantities at k^prime are set equal to quantities at k. k-dependent scalars and pointers associated phkxred=exp(.k.xred) for each atom [ham%gbound_k]=G sphere boundary, for each plane wave [ham%ph3d_k]=3-dim structure factors, for each atom and plane wave
SOURCE
990 subroutine load_k_hamiltonian(ham,ffnl_k,fockACE_k,gbound_k,istwf_k,kinpw_k,& 991 kg_k,kpg_k,kpt_k,npw_k,npw_fft_k,ph3d_k,& 992 compute_gbound,compute_ph3d) 993 994 !Arguments ------------------------------------ 995 !scalars 996 integer,intent(in),optional :: npw_k,npw_fft_k,istwf_k 997 logical,intent(in),optional :: compute_gbound,compute_ph3d 998 class(gs_hamiltonian_type),intent(inout),target :: ham 999 !arrays 1000 integer,intent(in),optional,target :: gbound_k(:,:),kg_k(:,:) 1001 real(dp),intent(in),optional :: kpt_k(3) 1002 real(dp),intent(in),optional,target :: ffnl_k(:,:,:,:),kinpw_k(:),kpg_k(:,:),ph3d_k(:,:,:) 1003 type(fock_ACE_type),intent(in),optional,target :: fockACE_k 1004 1005 !Local variables------------------------------- 1006 !scalars 1007 integer :: iat,iatom 1008 logical :: compute_gbound_ 1009 real(dp) :: arg 1010 !character(len=500) :: msg 1011 1012 ! ************************************************************************* 1013 1014 DBG_ENTER("COLL") 1015 1016 !@gs_hamiltonian_type 1017 1018 !k-dependent scalars 1019 if (present(kpt_k)) then 1020 ham%kpt_k(:) = kpt_k(:) 1021 ham%kpt_kp(:) = kpt_k(:) 1022 end if 1023 if (present(istwf_k)) then 1024 ham%istwf_k = istwf_k 1025 ham%istwf_kp = istwf_k 1026 end if 1027 if (present(npw_k)) then 1028 ham%npw_k = npw_k 1029 ham%npw_kp = npw_k 1030 end if 1031 if (present(npw_fft_k)) then 1032 ham%npw_fft_k = npw_fft_k 1033 ham%npw_fft_kp = npw_fft_k 1034 else if (present(npw_k)) then 1035 ham%npw_fft_k = npw_k 1036 ham%npw_fft_kp = npw_k 1037 end if 1038 1039 !Pointers to k-dependent quantitites 1040 if (present(kinpw_k)) then 1041 ham%kinpw_k => kinpw_k 1042 ham%kinpw_kp => kinpw_k 1043 end if 1044 if (present(kg_k)) then 1045 ham%kg_k => kg_k 1046 ham%kg_kp => kg_k 1047 end if 1048 if (present(kpg_k)) then 1049 ham%kpg_k => kpg_k 1050 ham%kpg_kp => kpg_k 1051 end if 1052 if (present(ffnl_k)) then 1053 ham%ffnl_k => ffnl_k 1054 ham%ffnl_kp => ffnl_k 1055 end if 1056 if (present(ph3d_k)) then 1057 ham%ph3d_k => ph3d_k 1058 ham%ph3d_kp => ph3d_k 1059 end if 1060 if (present(fockACE_k)) then 1061 ham%fockACE_k => fockACE_k 1062 end if 1063 !Compute exp(i.k.R) for each atom 1064 if (present(kpt_k)) then 1065 if (associated(Ham%phkpxred).and.(.not.associated(Ham%phkpxred,Ham%phkxred))) then 1066 ABI_FREE(Ham%phkpxred) 1067 end if 1068 ABI_SFREE(ham%phkxred) 1069 ABI_MALLOC(ham%phkxred,(2,ham%natom)) 1070 do iat=1,ham%natom 1071 iatom=ham%atindx(iat) 1072 arg=two_pi*DOT_PRODUCT(kpt_k,ham%xred(:,iat)) 1073 ham%phkxred(1,iatom)=DCOS(arg) 1074 ham%phkxred(2,iatom)=DSIN(arg) 1075 end do 1076 ham%phkpxred => ham%phkxred 1077 end if 1078 1079 !Compute or copy G sphere boundary at k+g 1080 compute_gbound_=.false.;if (present(compute_gbound)) compute_gbound_=compute_gbound 1081 if (present(gbound_k)) compute_gbound_=.true. 1082 if (compute_gbound_) then 1083 if (associated(Ham%gbound_kp,Ham%gbound_k)) then 1084 nullify(Ham%gbound_kp) 1085 else if (associated(Ham%gbound_kp)) then 1086 ABI_FREE(Ham%gbound_kp) 1087 end if 1088 ABI_SFREE(ham%gbound_k) 1089 end if 1090 if (.not.allocated(ham%gbound_k)) then 1091 ABI_MALLOC(ham%gbound_k,(2*ham%mgfft+8,2)) 1092 ham%gbound_k(:,:)=0 1093 ham%gbound_kp => ham%gbound_k 1094 end if 1095 if (compute_gbound_) then 1096 if (present(gbound_k)) then 1097 ham%gbound_k(:,:)=gbound_k(:,:) 1098 else 1099 if (.not.associated(ham%kg_k)) then 1100 ABI_BUG('Something is missing for gbound_k computation!') 1101 end if 1102 !write(std_out,*)"About to call sphereboundary" 1103 !write(std_out,*)"size(kg_k), npw_k, mgfft",size(ham%kg_k, dim=2), ham%npw_k, ham%mgfft 1104 call sphereboundary(ham%gbound_k,ham%istwf_k,ham%kg_k,ham%mgfft,ham%npw_k) 1105 end if 1106 ham%gbound_kp => ham%gbound_k 1107 end if 1108 1109 !Compute 3D structure factors for each atom at k+g 1110 if (present(compute_ph3d).and.present(ph3d_k)) then 1111 if (compute_ph3d.and.ham%nloalg(2)>0) then 1112 if ((.not.allocated(ham%phkxred)).or.(.not.associated(ham%kg_k)).or.& 1113 (.not.associated(ham%ph3d_k))) then 1114 ABI_BUG('Something is missing for ph3d_k computation!') 1115 end if 1116 call ph1d3d(1,ham%natom,ham%kg_k,ham%matblk,ham%natom,ham%npw_k,ham%ngfft(1),& 1117 ham%ngfft(2),ham%ngfft(3),ham%phkxred,ham%ph1d,ham%ph3d_k) 1118 end if 1119 end if 1120 1121 DBG_EXIT("COLL") 1122 1123 end subroutine load_k_hamiltonian
m_hamiltonian/load_k_rf_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
load_k_rf_hamiltonian
FUNCTION
Setup of the k-dependent part of the 1st- and 2nd- order Hamiltonian
INPUTS
[dkinpw_k]=1st derivative of the (modified) kinetic energy for each plane wave [ddkinpw_k]=2nd derivative of the (modified) kinetic energy for each plane wave [npw_k]=number of plane waves
SIDE EFFECTS
rf_Ham<rf_hamiltonian_type>=structured datatype completed with k-dependent quantitites. Quantities at k^prime are set equal to quantities at k.
SOURCE
1812 subroutine load_k_rf_hamiltonian(rf_Ham,dkinpw_k,ddkinpw_k,npw_k) 1813 1814 !Arguments ------------------------------------ 1815 !scalars 1816 integer,intent(in),optional :: npw_k 1817 class(rf_hamiltonian_type),intent(inout),target :: rf_Ham 1818 !arrays 1819 real(dp),intent(in),optional,target :: dkinpw_k(:),ddkinpw_k(:) 1820 1821 !Local variables------------------------------- 1822 1823 ! ************************************************************************* 1824 1825 DBG_ENTER("COLL") 1826 1827 !@gs_hamiltonian_type 1828 1829 !k-dependent scalars 1830 if (present(npw_k)) then 1831 rf_Ham%npw_k = npw_k 1832 rf_Ham%npw_kp = npw_k 1833 end if 1834 1835 !Pointers to k-dependent quantitites 1836 if (present(dkinpw_k)) then 1837 rf_Ham%dkinpw_k => dkinpw_k 1838 rf_Ham%dkinpw_kp => dkinpw_k 1839 end if 1840 if (present(ddkinpw_k)) then 1841 rf_Ham%ddkinpw_k => ddkinpw_k 1842 rf_Ham%ddkinpw_kp => ddkinpw_k 1843 end if 1844 1845 DBG_EXIT("COLL") 1846 1847 end subroutine load_k_rf_hamiltonian
m_hamiltonian/load_kprime_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
load_kprime_hamiltonian
FUNCTION
Setup of the k^prime-dependent part of the Hamiltonian H_k_k^prime
INPUTS
[compute_gbound]=flag. if true, G sphere boundary is computed here [compute_ph3d]=flag. if true, 3D structure factors are computed here (only if nloalg(2)>0) [gbound_kp]=G sphere boundary (not compatible with compute_gbound=TRUE) [ffnl_kp]=nonlocal form factors on basis sphere [istwf_kp]=parameter that describes the storage of wfs [kinpw_kp]=(modified) kinetic energy for each plane wave [kg_kp]=planewave reduced coordinates in basis sphere (g vectors) [kpg_kp]=(k+g) vectors in reciprocal space [kpt_kp]=k point coordinates [npw_kp]=number of plane waves (processed by current proc when band-FFT parallelism is on) [npw_fft_kp]=number of plane waves used to apply Hamiltonian (in the "FFT" configuration) [ph3d_kp]=3-dim structure factors, for each atom and plane wave
SIDE EFFECTS
ham<gs_hamiltonian_type>=structured datatype completed with k^prime-dependent quantitites. k^prime-dependent scalars and pointers associated phkpxred=exp(.k^prime.xred) for each atom [ham%gbound_kp]=G sphere boundary, for each plane wave [ham%ph3d_kp]=3-dim structure factors at k^prime at k, for each atom and plane wave
SOURCE
1158 subroutine load_kprime_hamiltonian(ham,ffnl_kp,gbound_kp,istwf_kp,kinpw_kp,& 1159 kg_kp,kpg_kp,kpt_kp,npw_kp,npw_fft_kp,& 1160 ph3d_kp,compute_gbound,compute_ph3d) 1161 1162 !Arguments ------------------------------------ 1163 !scalars 1164 integer,intent(in),optional :: npw_kp,npw_fft_kp,istwf_kp 1165 logical,intent(in),optional :: compute_gbound,compute_ph3d 1166 class(gs_hamiltonian_type),intent(inout),target :: ham 1167 !arrays 1168 integer,intent(in),optional,target :: gbound_kp(:,:),kg_kp(:,:) 1169 real(dp),intent(in),optional :: kpt_kp(3) 1170 real(dp),intent(in),optional,target :: ffnl_kp(:,:,:,:),kinpw_kp(:),kpg_kp(:,:),ph3d_kp(:,:,:) 1171 1172 !Local variables------------------------------- 1173 !scalars 1174 integer :: iat,iatom 1175 logical :: compute_gbound_ 1176 real(dp) :: arg 1177 !character(len=500) :: msg 1178 1179 ! ************************************************************************* 1180 1181 DBG_ENTER("COLL") 1182 1183 !@gs_hamiltonian_type 1184 1185 !k-dependent scalars 1186 if (present(kpt_kp)) ham%kpt_kp(:)= kpt_kp(:) 1187 if (present(istwf_kp)) ham%istwf_kp = istwf_kp 1188 if (present(npw_kp)) ham%npw_kp = npw_kp 1189 if (present(npw_fft_kp)) then 1190 ham%npw_fft_kp = npw_fft_kp 1191 else if (present(npw_kp)) then 1192 ham%npw_fft_kp = npw_kp 1193 end if 1194 1195 !Pointers to k-dependent quantitites 1196 if (present(kinpw_kp)) ham%kinpw_kp => kinpw_kp 1197 if (present(kg_kp)) ham%kg_kp => kg_kp 1198 if (present(kpg_kp)) ham%kpg_kp => kpg_kp 1199 if (present(ffnl_kp)) ham%ffnl_kp => ffnl_kp 1200 if (present(ph3d_kp)) ham%ph3d_kp => ph3d_kp 1201 1202 !Compute exp(i.k^prime.R) for each atom 1203 if (present(kpt_kp)) then 1204 if (associated(ham%phkpxred,ham%phkxred)) then 1205 nullify(ham%phkpxred) 1206 else if (associated(ham%phkpxred)) then 1207 ABI_FREE(ham%phkpxred) 1208 end if 1209 ABI_MALLOC(ham%phkpxred,(2,ham%natom)) 1210 do iat=1,ham%natom 1211 iatom=ham%atindx(iat) 1212 arg=two_pi*DOT_PRODUCT(kpt_kp,ham%xred(:,iat)) 1213 ham%phkpxred(1,iatom)=DCOS(arg) 1214 ham%phkpxred(2,iatom)=DSIN(arg) 1215 end do 1216 end if 1217 1218 !Compute or copy G sphere boundary at k^prime+g 1219 compute_gbound_=.false. 1220 if (present(kpt_kp).and.present(compute_gbound)) compute_gbound_=compute_gbound 1221 if (present(gbound_kp)) compute_gbound_=.true. 1222 if (compute_gbound_) then 1223 if (associated(ham%gbound_kp,ham%gbound_k)) then 1224 nullify(ham%gbound_kp) 1225 else if (associated(ham%gbound_kp)) then 1226 ABI_FREE(ham%gbound_kp) 1227 end if 1228 if (present(gbound_kp)) then 1229 ham%gbound_kp(:,:)=gbound_kp(:,:) 1230 else 1231 if (.not.associated(ham%kg_kp)) then 1232 ABI_BUG('Something is missing for gbound_kp computation!') 1233 end if 1234 ABI_MALLOC(ham%gbound_kp,(2*ham%mgfft+8,2)) 1235 call sphereboundary(ham%gbound_kp,ham%istwf_kp,ham%kg_kp,ham%mgfft,ham%npw_kp) 1236 end if 1237 end if 1238 1239 !Compute 3D structure factors for each atom at k^prime+g 1240 if (present(compute_ph3d).and.present(ph3d_kp)) then 1241 if (compute_ph3d.and.ham%nloalg(2)>0) then 1242 if ((.not.associated(ham%phkpxred)).or.(.not.associated(ham%kg_kp)).or.& 1243 & (.not.associated(ham%ph3d_kp))) then 1244 ABI_BUG('Something is missing for ph3d_kp computation!') 1245 end if 1246 call ph1d3d(1,ham%natom,ham%kg_kp,ham%matblk,ham%natom,ham%npw_kp,ham%ngfft(1),& 1247 & ham%ngfft(2),ham%ngfft(3),ham%phkpxred,ham%ph1d,ham%ph3d_kp) 1248 end if 1249 end if 1250 1251 DBG_EXIT("COLL") 1252 1253 end subroutine load_kprime_hamiltonian
m_hamiltonian/load_spin_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
load_spin_hamiltonian
INPUTS
isppol=index of current spin [vectornd(n4,n5,n6,nvloc,3)]=optional, vector potential of nuclear magnetic dipoles in real space [vlocal(n4,n5,n6,nvloc)]=optional, local potential in real space [vxctaulocal(n4,n5,n6,nvloc,4)]=optional, derivative of XC energy density with respect to kinetic energy density in real space [with_nonlocal]=optional, true if non-local factors have to be loaded
FUNCTION
Setup of the spin-dependent part of the GS Hamiltonian.
SIDE EFFECTS
Ham<gs_hamiltonian_type>=Structured datatype initialization phase: * Quantities that depend spin are initialized.
SOURCE
1456 subroutine load_spin_hamiltonian(Ham,isppol,vectornd,vlocal,vxctaulocal,with_nonlocal) 1457 1458 !Arguments ------------------------------------ 1459 !scalars 1460 integer,intent(in) :: isppol 1461 logical,optional,intent(in) :: with_nonlocal 1462 class(gs_hamiltonian_type),intent(inout),target :: Ham 1463 !arrays 1464 real(dp),optional,intent(in),target :: vectornd(:,:,:,:,:) 1465 real(dp),optional,intent(in),target :: vlocal(:,:,:,:),vxctaulocal(:,:,:,:,:) 1466 1467 !Local variables------------------------------- 1468 !scalars 1469 integer :: jsppol 1470 1471 ! ************************************************************************* 1472 1473 DBG_ENTER("COLL") 1474 1475 !@gs_hamiltonian_type 1476 if (present(vlocal)) then 1477 ABI_CHECK(size(vlocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc, "Wrong vlocal") 1478 Ham%vlocal => vlocal 1479 end if 1480 if (present(vxctaulocal)) then 1481 ABI_CHECK(size(vxctaulocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*4, "Wrong vxctaulocal") 1482 Ham%vxctaulocal => vxctaulocal 1483 end if 1484 if (present(vectornd)) then 1485 ABI_CHECK(size(vectornd)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*3, "Wrong vectornd") 1486 Ham%vectornd => vectornd 1487 end if 1488 1489 ! Retrieve non-local factors for this spin component 1490 if (present(with_nonlocal)) then 1491 if (with_nonlocal) then 1492 jsppol=min(isppol,size(Ham%ekb_spin,5)) 1493 if (jsppol>0) Ham%ekb => Ham%ekb_spin(:,:,:,:,jsppol) 1494 end if 1495 end if 1496 1497 ! Update enl and sij on GPU 1498 #if defined HAVE_GPU_CUDA 1499 if (Ham%gpu_option==ABI_GPU_LEGACY .or. Ham%gpu_option==ABI_GPU_KOKKOS) then 1500 call gpu_update_ham_data(& 1501 & Ham%ekb(:,:,:,1), INT(size(Ham%ekb), c_int64_t), & 1502 & Ham%sij, INT(size(Ham%sij), c_int64_t), & 1503 & Ham%gprimd, INT(size(Ham%gprimd),c_int64_t)) 1504 end if 1505 #endif 1506 1507 DBG_EXIT("COLL") 1508 1509 end subroutine load_spin_hamiltonian
m_hamiltonian/load_spin_rf_hamiltonian [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
load_spin_rf_hamiltonian
FUNCTION
Setup of the spin-dependent part of the 1st- and 2nd- order Hamiltonian.
INPUTS
isppol=index of current spin [vectornd(n4,n5,n6,nvloc)]=optional, vector potential of nuclear magnetic dipoles in real space in ddk direction idir [vlocal1(cplex*n4,n5,n6,nvloc)]=optional, 1st-order local potential in real space [vxctaulocal(n4,n5,n6,nvloc,4)]=optional, deriv of e_XC wrt kin energy, for mGGA [with_nonlocal]=optional, true if non-local factors have to be loaded
SIDE EFFECTS
rf_Ham<rf_hamiltonian_type>=Structured datatype initialization phase: * Quantities that depend on spin are initialized.
SOURCE
1736 subroutine load_spin_rf_hamiltonian(rf_Ham,isppol,vectornd,vlocal1,vxctaulocal,with_nonlocal) 1737 1738 !Arguments ------------------------------------ 1739 !scalars 1740 integer,intent(in) :: isppol 1741 logical,optional,intent(in) :: with_nonlocal 1742 class(rf_hamiltonian_type),intent(inout),target :: rf_Ham 1743 !arrays 1744 real(dp),optional,target,intent(in) :: vlocal1(:,:,:,:) 1745 real(dp),optional,target,intent(in) :: vectornd(:,:,:,:) 1746 real(dp),optional,target,intent(in) :: vxctaulocal(:,:,:,:,:) 1747 1748 !Local variables------------------------------- 1749 !scalars 1750 integer :: jsppol 1751 1752 ! ************************************************************************* 1753 1754 DBG_ENTER("COLL") 1755 1756 !@rf_hamiltonian_type 1757 1758 if (present(vlocal1)) then 1759 ABI_CHECK(size(vlocal1)==rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc,"Wrong vlocal1") 1760 rf_Ham%vlocal1 => vlocal1 1761 end if 1762 1763 if (present(vectornd)) then 1764 ABI_CHECK(size(vectornd)==rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc,"Wrong vectornd") 1765 rf_Ham%vectornd => vectornd 1766 end if 1767 1768 if (present(vxctaulocal)) then 1769 ABI_CHECK(size(vxctaulocal)==rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc*4,"Wrong vxctaulocal") 1770 rf_Ham%vxctaulocal => vxctaulocal 1771 end if 1772 1773 ! Retrieve non-local factors for this spin component 1774 if (present(with_nonlocal)) then 1775 if (with_nonlocal) then 1776 if (size(rf_Ham%e1kbfr_spin)>0) then 1777 jsppol=min(isppol,size(rf_Ham%e1kbfr_spin,5)) 1778 if (jsppol>0) rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsppol) 1779 end if 1780 if (size(rf_Ham%e1kbsc_spin)>0) then 1781 jsppol=min(isppol,size(rf_Ham%e1kbsc_spin,5)) 1782 if (jsppol>0) rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsppol) 1783 end if 1784 end if 1785 end if 1786 1787 DBG_EXIT("COLL") 1788 1789 end subroutine load_spin_rf_hamiltonian
m_hamiltonian/pawdij2e1kb [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
pawdij2e1kb
FUNCTION
Transfer PAW Dij (on-site RF Hamiltonian) values from paw_ij datastructure to e1kb array
INPUTS
OUTPUT
SOURCE
1944 subroutine pawdij2e1kb(paw_ij1,isppol,comm_atom,mpi_atmtab,e1kbfr,e1kbsc) 1945 1946 !Arguments ------------------------------------ 1947 !scalars 1948 integer,intent(in) :: isppol,comm_atom 1949 !arrays 1950 integer,intent(in),optional,target :: mpi_atmtab(:) 1951 real(dp),optional,intent(out) :: e1kbfr(:,:,:,:),e1kbsc(:,:,:,:) 1952 type(paw_ij_type),intent(in) :: paw_ij1(:) 1953 1954 !Local variables------------------------------- 1955 !scalars 1956 integer :: dimdij1,dime1kb1,dime1kb3,dime1kb4,iatom,iatom_tot,ierr,isp,ispden 1957 integer :: my_natom,natom,qphase 1958 logical :: my_atmtab_allocated,paral_atom 1959 !arrays 1960 integer,pointer :: my_atmtab(:) 1961 1962 ! ************************************************************************* 1963 1964 DBG_ENTER("COLL") 1965 1966 if ((.not.present(e1kbfr)).and.(.not.present(e1kbsc))) return 1967 if (present(e1kbfr)) then 1968 e1kbfr=zero ; natom=size(e1kbfr,2) 1969 end if 1970 if (present(e1kbsc)) then 1971 e1kbsc=zero ; natom=size(e1kbsc,2) 1972 end if 1973 1974 !Set up parallelism over atoms 1975 my_natom=size(paw_ij1) ; paral_atom=(xmpi_comm_size(comm_atom)>1) 1976 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 1977 call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 1978 1979 !Retrieve 1st-order PAW Dij coefficients for this spin component (frozen) 1980 if (my_natom>0.and.present(e1kbfr)) then 1981 if (allocated(paw_ij1(1)%dijfr)) then 1982 dime1kb1=size(e1kbfr,1) ; dime1kb3=size(e1kbfr,3) ; dime1kb4=size(e1kbfr,4) 1983 ABI_CHECK(paw_ij1(1)%qphase==dime1kb4,'BUG in pawdij2e1kb (1)!') 1984 do ispden=1,dime1kb3 1985 isp=isppol;if (dime1kb3==4) isp=ispden 1986 do iatom=1,my_natom 1987 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 1988 qphase=paw_ij1(iatom)%qphase 1989 dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size 1990 ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!') 1991 e1kbfr(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dijfr(1:dimdij1,isp) 1992 if (qphase==2) e1kbfr(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp) 1993 end do 1994 end do 1995 end if 1996 end if 1997 1998 !Retrieve 1st-order PAW Dij coefficients for this spin component (self-consistent) 1999 if (my_natom>0.and.present(e1kbsc)) then 2000 if (allocated(paw_ij1(1)%dijfr).and.allocated(paw_ij1(1)%dij)) then 2001 dime1kb1=size(e1kbsc,1) ; dime1kb3=size(e1kbsc,3) ; dime1kb4=size(e1kbsc,4) 2002 ABI_CHECK(paw_ij1(1)%qphase==dime1kb4,'BUG in pawdij2e1kb (1)!') 2003 do ispden=1,dime1kb3 2004 isp=isppol;if (dime1kb3==4) isp=ispden 2005 do iatom=1,my_natom 2006 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 2007 qphase=paw_ij1(iatom)%qphase 2008 dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size 2009 ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!') 2010 e1kbsc(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dij (1:dimdij1,isp) & 2011 & -paw_ij1(iatom)%dijfr(1:dimdij1,isp) 2012 if (qphase==2) e1kbsc(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dij (dimdij1+1:2*dimdij1,isp) & 2013 & -paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp) 2014 end do 2015 end do 2016 end if 2017 end if 2018 2019 ! Communication in case of distribution over atomic sites 2020 if (paral_atom) then 2021 if (present(e1kbfr)) call xmpi_sum(e1kbfr,comm_atom,ierr) 2022 if (present(e1kbsc)) call xmpi_sum(e1kbsc,comm_atom,ierr) 2023 end if 2024 2025 !Destroy atom table used for parallelism 2026 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 2027 2028 DBG_EXIT("COLL") 2029 2030 end subroutine pawdij2e1kb
m_hamiltonian/pawdij2ekb [ Functions ]
[ Top ] [ m_hamiltonian ] [ Functions ]
NAME
pawdij2ekb
FUNCTION
Transfer PAW Dij (on-site GS Hamiltonian) values from paw_ij datastructure to ekb array
INPUTS
OUTPUT
SOURCE
1866 subroutine pawdij2ekb(ekb,paw_ij,isppol,comm_atom,mpi_atmtab) 1867 1868 !Arguments ------------------------------------ 1869 !scalars 1870 integer,intent(in) :: isppol,comm_atom 1871 !arrays 1872 integer,intent(in),optional,target :: mpi_atmtab(:) 1873 real(dp),intent(out) :: ekb(:,:,:,:) 1874 type(paw_ij_type),intent(in) :: paw_ij(:) 1875 1876 !Local variables------------------------------- 1877 !scalars 1878 integer :: dimdij,dimekb1,dimekb3,dimekb4,iatom,iatom_tot,ierr,ii,isp,ispden,my_natom,natom,qphase 1879 logical :: my_atmtab_allocated,paral_atom 1880 !arrays 1881 integer,pointer :: my_atmtab(:) 1882 1883 ! ************************************************************************* 1884 1885 DBG_ENTER("COLL") 1886 1887 ekb=zero 1888 1889 !Set up parallelism over atoms 1890 natom=size(ekb,2); my_natom=size(paw_ij) 1891 paral_atom=(xmpi_comm_size(comm_atom)>1) 1892 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 1893 call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 1894 1895 !Retrieve PAW Dij coefficients for this spin component 1896 if (my_natom>0) then 1897 if (allocated(paw_ij(1)%dij)) then 1898 dimekb1=size(ekb,1) ; dimekb3=size(ekb,3) ; dimekb4=size(ekb,4) 1899 qphase=paw_ij(1)%qphase 1900 ABI_CHECK(qphase<=dimekb4,'paw_ij%qphase>dimekb4!') 1901 do ii=1,qphase 1902 do ispden=1,dimekb3 1903 isp=isppol; if (dimekb3==4) isp=ispden 1904 do iatom=1,my_natom 1905 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 1906 dimdij=paw_ij(iatom)%cplex_dij*paw_ij(iatom)%lmn2_size 1907 ABI_CHECK(dimdij<=dimekb1,'Size of paw_ij%dij>dimekb1!') 1908 ekb(1:dimdij,iatom_tot,ispden,ii)=paw_ij(iatom)%dij(1+(ii-1)*dimdij:ii*dimdij,isp) 1909 end do 1910 end do 1911 end do 1912 end if 1913 end if 1914 1915 !Communication in case of distribution over atomic sites 1916 if (paral_atom) then 1917 call xmpi_sum(ekb,comm_atom,ierr) 1918 end if 1919 1920 !Destroy atom table used for parallelism 1921 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 1922 1923 DBG_EXIT("COLL") 1924 1925 end subroutine pawdij2ekb
m_hamiltonian/rf_hamiltonian_type [ Types ]
[ Top ] [ m_hamiltonian ] [ Types ]
NAME
rf_hamiltonian_type
FUNCTION
This datastructure contains few data about one 1st-order Hamiltonian, needed in the "getgh1c" routine, that apply the 1st-order Hamiltonian on a wavefunction.
SOURCE
463 type,public :: rf_hamiltonian_type 464 465 ! ===== Integer scalars 466 467 integer :: cplex 468 ! if 1, real space 1-order functions on FFT grid are REAL; if 2, COMPLEX 469 470 integer :: dime1kb1 471 ! First dimension of E1kb, derivative of Ekb with respect to a perturbation 472 473 integer :: dime1kb2 474 ! Second dimension of E1kb, derivative of Ekb with respect to a perturbation 475 ! NCPP: dime1kb2=ntypat, PAW: dime1kb2=natom 476 477 integer :: npw_k 478 ! number of plane waves at k 479 480 integer :: npw_kp 481 ! number of plane waves at k^prime 482 483 integer:: nspinor 484 ! Number of spinorial components 485 486 integer :: nsppol 487 ! Total number of spin components (1=non-polarized, 2=polarized) 488 489 integer :: nvloc 490 ! Number of components of vloc 491 ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4 492 493 integer :: n4,n5,n6 494 ! same as ngfft(4:6) 495 496 ! ===== Real arrays 497 498 real(dp), allocatable :: e1kbfr_spin(:,:,:,:,:) 499 ! e1kbfr_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol) 500 ! Contains the values of e1kbfr array for all spins treated by current process 501 ! See e1kbfr description ; e1kbfr is pointer to e1kbfr_spin(:,:,:,:,isppol) 502 503 real(dp), allocatable :: e1kbsc_spin(:,:,:,:,:) 504 ! e1kbsc_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol) 505 ! Contains the values of e1kbsc array for all spins treated by current process 506 ! See e1kbsc description ; e1kbsc is pointer to e1kbsc_spin(:,:,:,:,isppol) 507 508 ! ===== Real pointers 509 510 real(dp), pointer :: dkinpw_k(:) => null() 511 ! dkinpw_k(npw_k) 512 ! 1st derivative of the (modified) kinetic energy for each plane wave at k 513 514 real(dp), pointer :: dkinpw_kp(:) => null() 515 ! dkinpw_kp(npw_kp) 516 ! 1st derivative of the (modified) kinetic energy for each plane wave at k^prime 517 518 real(dp), pointer :: ddkinpw_k(:) => null() 519 ! ddkinpw_k(npw_k) 520 ! 2nd derivative of the (modified) kinetic energy for each plane wave at k 521 522 real(dp), pointer :: ddkinpw_kp(:) => null() 523 ! ddkinpw_kp(npw_kp) 524 ! 2nd derivative of the (modified) kinetic energy for each plane wave at k^prime 525 526 real(dp), pointer :: e1kbfr(:,:,:,:) => null() 527 ! Frozen part of 1st derivative of ekb for the considered perturbation 528 ! (part not depending on VHxc^(1)) 529 ! e1kbfr(dime1kb1,dime1kb2,nspinor**2,cplex) 530 ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol) 531 532 real(dp), ABI_CONTIGUOUS pointer :: e1kbsc(:,:,:,:) => null() 533 ! Self-consistent 1st derivative of ekb for the considered perturbation 534 ! (part depending only on self-consistent VHxc^(1)) 535 ! e1kbsc(dime1kb1,dime1kb2,nspinor**2,cplex) 536 ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol) 537 538 real(dp), pointer :: vectornd(:,:,:,:) => null() 539 ! vectornd(n4,n5,n6,nvloc) 540 ! vector potential of nuclear magnetic dipoles 541 ! in real space, on the augmented fft grid, in direction idir 542 ! (the ddk pert direction) 543 544 real(dp), pointer :: vlocal1(:,:,:,:) => null() 545 ! vlocal1(cplex*n4,n5,n6,nvloc) 546 ! 1st-order local potential in real space, on the augmented fft grid 547 548 real(dp), pointer :: vxctaulocal(:,:,:,:,:) => null() 549 ! vxctaulocal(n4,n5,n6,nvloc,4) 550 ! derivative of XC energy density with respect to kinetic energy density, 551 ! in real space, on the augmented fft grid 552 553 contains 554 procedure :: free => destroy_rf_hamiltonian 555 ! Free the memory in the RF Hamiltonian 556 557 procedure :: load_spin => load_spin_rf_hamiltonian 558 ! Setup of the spin-dependent part of the RF Hamiltonian. 559 560 procedure :: load_k => load_k_rf_hamiltonian 561 ! Setup of the k-dependent part of the RF Hamiltonian 562 563 end type rf_hamiltonian_type 564 565 public :: init_rf_hamiltonian ! Initialize the RF Hamiltonian