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
1977 subroutine gspot_transgrid_and_pack(isppol, usepaw, paral_kgb, nfft, ngfft, nfftf, & 1978 nspden, nvloc, ncomp, pawfgr, mpi_enreg, vtrial, vlocal) 1979 1980 !Arguments ------------------------------- 1981 integer,intent(in) :: isppol, nspden, ncomp, usepaw, paral_kgb, nfft, nfftf, nvloc 1982 type(pawfgr_type), intent(in) :: pawfgr 1983 type(MPI_type), intent(in) :: mpi_enreg 1984 !arrays 1985 integer,intent(in) :: ngfft(18) 1986 real(dp),intent(inout) :: vtrial(nfftf, nspden, ncomp) 1987 real(dp),intent(out) :: vlocal(ngfft(4), ngfft(5), ngfft(6), nvloc, ncomp) 1988 1989 1990 !Local variables------------------------------- 1991 !scalars 1992 integer :: n1,n2,n3,n4,n5,n6,ispden,ic 1993 real(dp) :: rhodum(1) 1994 real(dp),allocatable :: cgrvtrial(:,:), vlocal_tmp(:,:,:) 1995 1996 ! ************************************************************************* 1997 1998 ! Coarse mesh. 1999 n1=ngfft(1); n2=ngfft(2); n3=ngfft(3) 2000 n4=ngfft(4); n5=ngfft(5); n6=ngfft(6) 2001 2002 ! Set up local potential vlocal with proper dimensioning, from vtrial 2003 ! Also take into account the spin. 2004 if (nspden /= 4) then 2005 if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then 2006 ! Fine mesh == Coarse mesh. 2007 do ic=1,ncomp 2008 call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal(:,:,:,:,ic),2) 2009 end do 2010 else 2011 ! Transfer from fine mesh to coarse and then pack data 2012 ABI_MALLOC(cgrvtrial,(nfft,nspden)) 2013 do ic=1,ncomp 2014 call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,& 2015 rhodum,rhodum,cgrvtrial,vtrial(:,:,ic)) 2016 call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,& 2017 cgrvtrial,vlocal(:,:,:,:,ic),2) 2018 end do 2019 ABI_FREE(cgrvtrial) 2020 end if 2021 else 2022 ! nspden == 4. replace isppol by loop over ispden. 2023 ABI_MALLOC(vlocal_tmp, (n4,n5,n6)) 2024 if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then 2025 ! Fine mesh == Coarse mesh. 2026 do ic=1,ncomp 2027 do ispden=1,nspden 2028 call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal_tmp,2) 2029 vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:) 2030 end do 2031 end do 2032 else 2033 ! Transfer from fine mesh to coarse and then pack data 2034 ABI_MALLOC(cgrvtrial,(nfft,nspden)) 2035 do ic=1,ncomp 2036 call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,rhodum,rhodum,cgrvtrial,vtrial(:,:,ic)) 2037 do ispden=1,nspden 2038 call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,cgrvtrial,vlocal_tmp,2) 2039 vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:) 2040 end do 2041 end do 2042 ABI_FREE(cgrvtrial) 2043 end if 2044 ABI_FREE(vlocal_tmp) 2045 end if ! nspden 2046 2047 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-2022 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 defs_basis 33 use m_abicore 34 use m_errors 35 use m_xmpi 36 37 use defs_datatypes, only : pseudopotential_type 38 use defs_abitypes, only : MPI_type 39 use m_copy, only : addr_copy 40 use m_geometry, only : metric 41 use m_pawtab, only : pawtab_type 42 use m_pawfgr, only : pawfgr_type 43 use m_fftcore, only : sphereboundary 44 use m_fft, only : fftpac 45 use m_fourier_interpol, only : transgrid 46 use m_pawcprj, only : pawcprj_getdim 47 use m_paw_ij, only : paw_ij_type 48 use m_paral_atom, only : get_my_atmtab, free_my_atmtab 49 use m_electronpositron, only : electronpositron_type, electronpositron_calctype 50 use m_kg, only : ph1d3d, getph 51 use m_fock, only : fock_common_type, fock_BZ_type, fock_ACE_type, fock_type 52 53 #if defined HAVE_GPU_CUDA 54 use m_manage_cuda 55 #endif 56 57 #if defined HAVE_FC_ISO_C_BINDING 58 use, intrinsic :: iso_c_binding, only : c_ptr,c_loc,c_f_pointer,c_int32_t 59 #endif 60 61 implicit none 62 63 private 64 65 public :: pawdij2ekb 66 public :: pawdij2e1kb 67 public :: gspot_transgrid_and_pack ! Set up local potential vlocal on the coarse FFT mesh from vtrial on the fine mesh. 68 69 ! These constants select how H is applied in reciprocal space 70 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
1210 subroutine copy_hamiltonian(gs_hamk_out,gs_hamk_in) 1211 1212 !Arguments ------------------------------------ 1213 type(gs_hamiltonian_type),intent(in),target :: gs_hamk_in 1214 type(gs_hamiltonian_type),intent(out),target :: gs_hamk_out 1215 1216 !Local variables------------------------------- 1217 integer :: tmp2i(5) 1218 #if defined HAVE_FC_ISO_C_BINDING 1219 type(C_PTR) :: ham_ptr 1220 #endif 1221 1222 ! ************************************************************************* 1223 1224 DBG_ENTER("COLL") 1225 1226 !@gs_hamiltonian_type 1227 1228 gs_hamk_out%dimekb1 = gs_hamk_in%dimekb1 1229 gs_hamk_out%dimekb2 = gs_hamk_in%dimekb2 1230 gs_hamk_out%dimekbq = gs_hamk_in%dimekbq 1231 gs_hamk_out%istwf_k = gs_hamk_in%istwf_k 1232 gs_hamk_out%istwf_kp = gs_hamk_in%istwf_kp 1233 gs_hamk_out%lmnmax = gs_hamk_in%lmnmax 1234 gs_hamk_out%matblk = gs_hamk_in%matblk 1235 gs_hamk_out%mgfft = gs_hamk_in%mgfft 1236 gs_hamk_out%mpsang = gs_hamk_in%mpsang 1237 gs_hamk_out%mpssoang = gs_hamk_in%mpssoang 1238 gs_hamk_out%natom = gs_hamk_in%natom 1239 gs_hamk_out%nfft = gs_hamk_in%nfft 1240 gs_hamk_out%npw_k = gs_hamk_in%npw_k 1241 gs_hamk_out%npw_kp = gs_hamk_in%npw_kp 1242 gs_hamk_out%npw_fft_k = gs_hamk_in%npw_fft_k 1243 gs_hamk_out%npw_fft_kp = gs_hamk_in%npw_fft_kp 1244 gs_hamk_out%nspinor = gs_hamk_in%nspinor 1245 gs_hamk_out%nsppol = gs_hamk_in%nsppol 1246 gs_hamk_out%ntypat = gs_hamk_in%ntypat 1247 gs_hamk_out%nvloc = gs_hamk_in%nvloc 1248 gs_hamk_out%n4 = gs_hamk_in%n4 1249 gs_hamk_out%n5 = gs_hamk_in%n5 1250 gs_hamk_out%n6 = gs_hamk_in%n6 1251 gs_hamk_out%use_gpu_cuda = gs_hamk_in%use_gpu_cuda 1252 gs_hamk_out%usecprj = gs_hamk_in%usecprj 1253 gs_hamk_out%usepaw = gs_hamk_in%usepaw 1254 gs_hamk_out%useylm = gs_hamk_in%useylm 1255 gs_hamk_out%ngfft = gs_hamk_in%ngfft 1256 gs_hamk_out%nloalg = gs_hamk_in%nloalg 1257 gs_hamk_out%ucvol = gs_hamk_in%ucvol 1258 gs_hamk_out%gmet = gs_hamk_in%gmet 1259 gs_hamk_out%gprimd = gs_hamk_in%gprimd 1260 gs_hamk_out%kpt_k = gs_hamk_in%kpt_k 1261 gs_hamk_out%kpt_kp = gs_hamk_in%kpt_kp 1262 1263 ABI_MALLOC(gs_hamk_out%atindx,(gs_hamk_out%natom)) 1264 gs_hamk_out%atindx = gs_hamk_in%atindx 1265 ABI_MALLOC(gs_hamk_out%atindx1,(gs_hamk_out%natom)) 1266 gs_hamk_out%atindx1 = gs_hamk_in%atindx1 1267 ABI_MALLOC(gs_hamk_out%dimcprj,(gs_hamk_out%natom*gs_hamk_out%usepaw)) 1268 if (gs_hamk_out%usepaw==1) gs_hamk_out%dimcprj = gs_hamk_in%dimcprj 1269 ABI_MALLOC(gs_hamk_out%typat,(gs_hamk_out%natom)) 1270 gs_hamk_out%typat = gs_hamk_in%typat 1271 ABI_MALLOC(gs_hamk_out%gbound_k,(2*gs_hamk_out%mgfft+8,2)) 1272 gs_hamk_out%gbound_k = gs_hamk_in%gbound_k 1273 ABI_MALLOC(gs_hamk_out%indlmn,(6,gs_hamk_out%lmnmax,gs_hamk_out%ntypat)) 1274 gs_hamk_out%indlmn = gs_hamk_in%indlmn 1275 ABI_MALLOC(gs_hamk_out%nattyp,(gs_hamk_out%ntypat)) 1276 gs_hamk_out%nattyp = gs_hamk_in%nattyp 1277 ABI_MALLOC(gs_hamk_out%nucdipmom,(3,gs_hamk_out%natom)) 1278 gs_hamk_out%nucdipmom = gs_hamk_in%nucdipmom 1279 ABI_MALLOC(gs_hamk_out%phkxred,(2,gs_hamk_out%natom)) 1280 gs_hamk_out%phkxred = gs_hamk_in%phkxred 1281 ABI_MALLOC(gs_hamk_out%ph1d,(2,3*(2*gs_hamk_out%mgfft+1)*gs_hamk_out%natom)) 1282 gs_hamk_out%ph1d = gs_hamk_in%ph1d 1283 ABI_MALLOC(gs_hamk_out%pspso,(gs_hamk_out%ntypat)) 1284 gs_hamk_out%pspso = gs_hamk_in%pspso 1285 tmp2i(1:5)=shape(gs_hamk_in%ekb_spin) 1286 ABI_MALLOC(gs_hamk_out%ekb_spin,(tmp2i(1),tmp2i(2),tmp2i(3),tmp2i(4),tmp2i(5))) 1287 gs_hamk_out%ekb_spin = gs_hamk_in%ekb_spin 1288 gs_hamk_out%ekb => gs_hamk_out%ekb_spin(:,:,:,:,1) 1289 tmp2i(1:2)=shape(gs_hamk_in%sij) 1290 ABI_MALLOC(gs_hamk_out%sij,(tmp2i(1),tmp2i(2))) 1291 gs_hamk_out%sij = gs_hamk_in%sij 1292 1293 if (associated(gs_hamk_in%gbound_kp,gs_hamk_in%gbound_k)) then 1294 gs_hamk_out%gbound_kp => gs_hamk_out%gbound_k 1295 else 1296 ABI_MALLOC(gs_hamk_out%gbound_kp,(2,gs_hamk_out%natom)) 1297 gs_hamk_out%gbound_kp = gs_hamk_in%gbound_kp 1298 end if 1299 if (associated(gs_hamk_in%phkpxred,gs_hamk_in%phkxred)) then 1300 gs_hamk_out%phkpxred => gs_hamk_out%phkxred 1301 else 1302 ABI_MALLOC(gs_hamk_out%phkpxred,(2,gs_hamk_out%natom)) 1303 gs_hamk_out%phkpxred = gs_hamk_in%phkpxred 1304 end if 1305 1306 call addr_copy(gs_hamk_in%xred,gs_hamk_out%xred) 1307 call addr_copy(gs_hamk_in%vectornd,gs_hamk_out%vectornd) 1308 call addr_copy(gs_hamk_in%vlocal,gs_hamk_out%vlocal) 1309 call addr_copy(gs_hamk_in%vxctaulocal,gs_hamk_out%vxctaulocal) 1310 call addr_copy(gs_hamk_in%kinpw_k,gs_hamk_out%kinpw_k) 1311 call addr_copy(gs_hamk_in%kinpw_kp,gs_hamk_out%kinpw_kp) 1312 call addr_copy(gs_hamk_in%kg_k,gs_hamk_out%kg_k) 1313 call addr_copy(gs_hamk_in%kg_kp,gs_hamk_out%kg_kp) 1314 call addr_copy(gs_hamk_in%kpg_k,gs_hamk_out%kpg_k) 1315 call addr_copy(gs_hamk_in%kpg_kp,gs_hamk_out%kpg_kp) 1316 call addr_copy(gs_hamk_in%ffnl_k,gs_hamk_out%ffnl_k) 1317 call addr_copy(gs_hamk_in%ffnl_kp,gs_hamk_out%ffnl_kp) 1318 call addr_copy(gs_hamk_in%ph3d_k,gs_hamk_out%ph3d_k) 1319 call addr_copy(gs_hamk_in%ph3d_kp,gs_hamk_out%ph3d_kp) 1320 1321 !For pointers to structured datatypes, have to copy the address 1322 !manually because there is no generic addr_copy function for that 1323 if (associated(gs_hamk_in%fockcommon)) then 1324 #if defined HAVE_FC_ISO_C_BINDING 1325 ham_ptr=c_loc(gs_hamk_in%fockcommon) 1326 call c_f_pointer(ham_ptr,gs_hamk_out%fockcommon) 1327 #else 1328 gs_hamk_out%fockcommon=transfer(gs_hamk_in%fockcommon,gs_hamk_out%fockcommon) 1329 #endif 1330 else 1331 nullify(gs_hamk_out%fockcommon) 1332 end if 1333 if (associated(gs_hamk_in%fockbz)) then 1334 #if defined HAVE_FC_ISO_C_BINDING 1335 ham_ptr=c_loc(gs_hamk_in%fockbz) 1336 call c_f_pointer(ham_ptr,gs_hamk_out%fockbz) 1337 #else 1338 gs_hamk_out%fockbz=transfer(gs_hamk_in%fockbz,gs_hamk_out%fockbz) 1339 #endif 1340 else 1341 nullify(gs_hamk_out%fockbz) 1342 end if 1343 if (associated(gs_hamk_in%fockACE_k)) then 1344 #if defined HAVE_FC_ISO_C_BINDING 1345 ham_ptr=c_loc(gs_hamk_in%fockACE_k) 1346 call c_f_pointer(ham_ptr,gs_hamk_out%fockACE_k) 1347 #else 1348 gs_hamk_out%fockACE_k=transfer(gs_hamk_in%fockACE_k,gs_hamk_out%fockACE_k) 1349 #endif 1350 else 1351 nullify(gs_hamk_out%fockACE_k) 1352 end if 1353 1354 DBG_EXIT("COLL") 1355 1356 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
543 subroutine destroy_hamiltonian(Ham) 544 545 !Arguments ------------------------------------ 546 !scalars 547 class(gs_hamiltonian_type),intent(inout),target :: Ham 548 549 ! ************************************************************************* 550 551 DBG_ENTER("COLL") 552 553 !@gs_hamiltonian_type 554 555 ! Integer Pointers 556 if (associated(Ham%gbound_kp,Ham%gbound_k)) then 557 nullify(Ham%gbound_kp) 558 else if (associated(Ham%gbound_kp)) then 559 ABI_FREE(Ham%gbound_kp) 560 end if 561 562 ! Integer arrays 563 ABI_SFREE(Ham%atindx) 564 ABI_SFREE(Ham%atindx1) 565 ABI_SFREE(Ham%dimcprj) 566 ABI_SFREE(Ham%gbound_k) 567 ABI_SFREE(Ham%indlmn) 568 ABI_SFREE(Ham%nattyp) 569 ABI_SFREE(Ham%pspso) 570 ABI_SFREE(Ham%typat) 571 572 ! Real Pointers 573 if (associated(Ham%phkpxred,Ham%phkxred)) then 574 nullify(Ham%phkpxred) 575 else if (associated(Ham%phkpxred)) then 576 ABI_FREE(Ham%phkpxred) 577 end if 578 ABI_SFREE(Ham%phkxred) 579 if (associated(Ham%ekb)) nullify(Ham%ekb) 580 if (associated(Ham%vectornd)) nullify(Ham%vectornd) 581 if (associated(Ham%vlocal)) nullify(Ham%vlocal) 582 if (associated(Ham%vxctaulocal)) nullify(Ham%vxctaulocal) 583 if (associated(Ham%xred)) nullify(Ham%xred) 584 if (associated(Ham%kinpw_k)) nullify(Ham%kinpw_k) 585 if (associated(Ham%kinpw_kp)) nullify(Ham%kinpw_kp) 586 if (associated(Ham%kg_k)) nullify(Ham%kg_k) 587 if (associated(Ham%kg_kp)) nullify(Ham%kg_kp) 588 if (associated(Ham%kpg_k)) nullify(Ham%kpg_k) 589 if (associated(Ham%kpg_kp)) nullify(Ham%kpg_kp) 590 if (associated(Ham%ffnl_k)) nullify(Ham%ffnl_k) 591 if (associated(Ham%ffnl_kp)) nullify(Ham%ffnl_kp) 592 if (associated(Ham%ph3d_k)) nullify(Ham%ph3d_k) 593 if (associated(Ham%ph3d_kp)) nullify(Ham%ph3d_kp) 594 595 ! Real arrays 596 ABI_SFREE(Ham%ekb_spin) 597 ABI_SFREE(Ham%sij) 598 ABI_SFREE(Ham%nucdipmom) 599 ABI_SFREE(Ham%ph1d) 600 601 ! Structured datatype pointers 602 if (associated(Ham%fockcommon)) nullify(Ham%fockcommon) 603 if (associated(Ham%fockACE_k)) nullify(Ham%fockACE_k) 604 if (associated(Ham%fockbz)) nullify(Ham%fockbz) 605 #if defined HAVE_GPU_CUDA 606 if(Ham%use_gpu_cuda==1) then 607 call gpu_finalize_ham_data() 608 end if 609 #endif 610 611 DBG_EXIT("COLL") 612 613 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
1449 subroutine destroy_rf_hamiltonian(rf_Ham) 1450 1451 !Arguments ------------------------------------ 1452 !scalars 1453 class(rf_hamiltonian_type),intent(inout) :: rf_Ham 1454 1455 ! ************************************************************************* 1456 1457 DBG_ENTER("COLL") 1458 1459 !@rf_hamiltonian_type 1460 1461 ! Real arrays 1462 ABI_SFREE(rf_Ham%e1kbfr_spin) 1463 ABI_SFREE(rf_Ham%e1kbsc_spin) 1464 1465 ! Real pointers 1466 if (associated(rf_Ham%dkinpw_k)) nullify(rf_Ham%dkinpw_k) 1467 if (associated(rf_Ham%dkinpw_kp)) nullify(rf_Ham%dkinpw_kp) 1468 if (associated(rf_Ham%ddkinpw_k)) nullify(rf_Ham%ddkinpw_k) 1469 if (associated(rf_Ham%ddkinpw_kp)) nullify(rf_Ham%ddkinpw_kp) 1470 if (associated(rf_Ham%vectornd)) nullify(rf_Ham%vectornd) 1471 if (associated(rf_Ham%vlocal1)) nullify(rf_Ham%vlocal1) 1472 if (associated(rf_Ham%e1kbfr)) nullify(rf_Ham%e1kbfr) 1473 if (associated(rf_Ham%e1kbsc)) nullify(rf_Ham%e1kbsc) 1474 1475 DBG_EXIT("COLL") 1476 1477 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
90 type,public :: gs_hamiltonian_type 91 92 ! ===== Integer scalars 93 94 integer :: dimekb1 95 ! First dimension of Ekb 96 ! Same as psps%dimekb 97 ! ->Norm conserving : Max. number of Kleinman-Bylander energies 98 ! for each atom type 99 ! dimekb1=lnmax 100 ! ->PAW : Max. number of Dij coefficients connecting projectors 101 ! for each atom 102 ! dimekb1=cplex_dij*lmnmax*(lmnmax+1)/2 103 104 integer :: dimekb2 105 ! Second dimension of Ekb 106 ! ->Norm conserving psps: dimekb2=ntypat 107 ! ->PAW : dimekb2=natom 108 109 integer :: dimekbq 110 ! Fourth dimension of Ekb 111 ! 2 if Ekb factors contain a exp(-iqR) phase, 1 otherwise 112 113 integer :: istwf_k 114 ! option parameter that describes the storage of wfs at k 115 116 integer :: istwf_kp 117 ! option parameter that describes the storage of wfs at k^prime 118 119 integer :: lmnmax 120 ! Maximum number of different l,m,n components over all types of psps. 121 ! same as dtset%lmnmax 122 123 integer :: matblk 124 ! dimension of the array ph3d 125 126 integer :: mgfft 127 ! maximum size for 1D FFTs (same as dtset%mgfft) 128 129 integer :: mpsang 130 ! Highest angular momentum of non-local projectors over all type of psps. 131 ! shifted by 1 : for all local psps, mpsang=0; for largest s, mpsang=1, 132 ! for largest p, mpsang=2; for largest d, mpsang=3; for largest f, mpsang=4 133 ! This gives also the number of non-local "channels" 134 ! same as psps%mpsang 135 136 integer :: mpssoang 137 ! Maximum number of channels, including those for treating the spin-orbit coupling 138 ! For NC pseudopotentials only: 139 ! when mpspso=1, mpssoang=mpsang 140 ! when mpspso=2, mpssoang=2*mpsang-1 141 ! For PAW: same as mpsang 142 ! same as psps%mpssoang 143 144 integer :: natom 145 ! The number of atoms for this dataset; same as dtset%natom 146 147 integer :: nfft 148 ! number of FFT grid points same as dtset%nfft 149 150 integer :: npw_k 151 ! number of plane waves at k 152 ! In case of band-FFT parallelism, npw_k is the number of plane waves 153 ! processed by current proc 154 155 integer :: npw_fft_k 156 ! number of plane waves at k used to apply Hamiltonian when band-FFT 157 ! parallelism is activated (i.e. data are distributed in the "FFT" configuration) 158 159 integer :: npw_kp 160 ! number of plane waves at k^prime 161 ! In case of band-FFT parallelism, npw_kp is the number of plane waves 162 ! processed by current proc 163 164 integer :: npw_fft_kp 165 ! number of plane waves at k^prime used to apply Hamiltonian when band-FFT 166 ! parallelism is activated (i.e. data are distributed in the "FFT" configuration) 167 168 integer :: nspinor 169 ! Number of spinorial components 170 171 integer :: nsppol 172 ! Total number of spin components (1=non-polarized, 2=polarized) 173 174 integer :: ntypat 175 ! Number of types of pseudopotentials same as dtset%ntypat 176 177 integer :: nvloc 178 ! Number of components of vloc 179 ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4 180 181 integer :: n4,n5,n6 182 ! same as ngfft(4:6) 183 184 integer :: use_gpu_cuda 185 ! governs wheter we do the hamiltonian calculation on gpu or not 186 187 integer :: usecprj 188 ! usecprj= 1 if cprj projected WF are stored in memory 189 ! = 0 if they are to be computed on the fly 190 191 integer :: usepaw 192 ! if usepaw=0 , use norm-conserving psps part of the code 193 ! is usepaw=1 , use paw part of the code 194 195 integer :: useylm 196 ! governs the way the nonlocal operator is to be applied: 197 ! 1=using Ylm, 0=using Legendre polynomials 198 199 ! ===== Integer arrays 200 201 integer, allocatable :: atindx(:) 202 ! atindx(natom) 203 ! index table for atoms (see gstate.f) 204 205 integer, allocatable :: atindx1(:) 206 ! atindx1(natom) 207 ! index table for atoms, inverse of atindx (see gstate.f) 208 209 integer, allocatable :: dimcprj(:) 210 ! dimcprj(natom*usepaw)=dimensions of array cprj 211 ! dimcprj(ia)=cprj(ia,:)%nlmn 212 ! atoms are ordered by atom-type 213 214 integer, allocatable :: gbound_k(:,:) 215 ! gbound_k(2*mgfft+8,2) 216 ! G sphere boundary, for each plane wave at k 217 218 integer(kind=c_int32_t), allocatable :: indlmn(:,:,:) 219 ! indlmn(6,lmnmax,ntypat) 220 ! For each type of psp, 221 ! array giving l,m,n,lm,ln,spin for i=ln (if useylm=0) 222 ! or i=lmn (if useylm=1) 223 224 integer, allocatable :: nattyp(:) 225 ! nattyp(ntypat) 226 ! # of atoms of each type 227 228 integer :: ngfft(18) 229 ! ngfft(1:3)=integer fft box dimensions 230 ! ngfft(4:6)=integer fft box dimensions, might be augmented for CPU speed 231 ! ngfft(7)=fftalg 232 ! ngfft(8)=fftalg 233 234 integer :: nloalg(3) 235 ! governs the choice of the algorithm for non-local operator same as dtset%nloalg 236 237 integer, allocatable :: pspso(:) 238 ! pspso(ntypat) 239 ! For each type of psp, 1 if no spin-orbit component is taken 240 ! into account, 2 if a spin-orbit component is used 241 ! Revelant for NC-psps and PAW 242 243 integer, allocatable :: typat(:) 244 ! typat(natom) 245 ! type of each atom 246 247 ! integer, allocatable :: indpw_k(:,:) 248 ! indpw_k(4,npw_fft_k) 249 ! array which gives fft box index for given basis sphere 250 251 ! Integer pointers 252 253 integer, ABI_CONTIGUOUS pointer :: gbound_kp(:,:) => null() 254 ! gbound_kp(2*mgfft+8,2) 255 ! G sphere boundary, for each plane wave at k^prime 256 257 integer, pointer :: kg_k(:,:) => null() 258 ! kg_k(3,npw_fft_k) 259 ! G vector coordinates with respect to reciprocal lattice translations 260 ! at k 261 262 integer, pointer :: kg_kp(:,:) => null() 263 ! kg_kp(3,npw_fft_kp) 264 ! G vector coordinates with respect to reciprocal lattice translations 265 ! at k^prime 266 267 ! ===== Real scalars 268 269 real(dp) :: ucvol 270 ! unit cell volume (Bohr**3) 271 272 ! ===== Real arrays 273 274 real(dp), allocatable :: ekb_spin(:,:,:,:,:) 275 ! ekb_spin(dimekb1,dimekb2,nspinor**2,dimekbq,my_nsppol) 276 ! Contains the values of ekb array for all spins treated by current process 277 ! See ekb description ; ekb is pointer to ekb_spin(:,:,:,:,my_isppol) 278 279 real(dp), allocatable :: sij(:,:) 280 ! sij(dimekb1,ntypat*usepaw) = overlap matrix for paw calculation 281 282 real(dp) :: gmet(3,3) 283 ! reciprocal space metric tensor in Bohr**-2 284 285 real(dp) :: gprimd(3,3) 286 ! dimensional reciprocal space primitive translations (Bohr^-1) 287 288 real(dp) :: kpt_k(3) 289 ! dimensionless k point coordinates wrt reciprocal lattice vectors 290 291 real(dp) :: kpt_kp(3) 292 ! dimensionless k^prime point coordinates wrt reciprocal lattice vectors 293 294 real(dp), allocatable :: nucdipmom(:,:) 295 ! nucdipmom(3,natom) 296 ! nuclear dipole moments at each atomic position 297 298 real(dp), allocatable :: ph1d(:,:) 299 ! ph1d(2,3*(2*mgfft+1)*natom) 300 ! 1-dim phase arrays for structure factor (see getph.f). 301 302 real(dp), allocatable :: phkxred(:,:) 303 ! phkxred(2,natom) 304 ! phase factors exp(2 pi k.xred) at k 305 306 ! ===== Real pointers 307 308 real(dp), ABI_CONTIGUOUS pointer :: ekb(:,:,:,:) => null() 309 ! ekb(dimekb1,dimekb2,nspinor**2,dimekbq) 310 ! ->Norm conserving : (Real) Kleinman-Bylander energies (hartree) 311 ! for number of basis functions (l,n) (lnmax) 312 ! and number of atom types (ntypat) 313 ! dimekb1=lnmax ; dimekb2=ntypat ; dimekbq=1 314 ! ->PAW : (Real, symmetric) Frozen part of Dij coefficients 315 ! to connect projectors 316 ! for number of basis functions (l,m,n) (lmnmax) 317 ! and number of atom (natom) 318 ! dimekb1=lmnmax*(lmnmax+1)/2 ; dimekb2=natom ; dimekbq=1 319 ! ekb is spin dependent in the case of PAW calculations. 320 ! For each spin component, ekb points to ekb_spin(:,:,:,:,my_isppol) 321 ! dimekbq=2 if Ekb factors contain a exp(-iqR) phase, dimekbq=1 otherwise 322 ! About the non-local factors symmetry: 323 ! - The lower triangular part of the Dij matrix can be deduced from the upper one 324 ! with the following relation: D^s2s1_ji = (D^s1s2_ij)^* 325 ! where s1,s2 are spinor components 326 327 real(dp), pointer :: ffnl_k(:,:,:,:) => null() 328 ! ffnl_k(npw_fft_k,2,dimffnl_k,ntypat) 329 ! nonlocal form factors at k 330 331 real(dp), pointer :: ffnl_kp(:,:,:,:) => null() 332 ! ffnl_kp(npw_fft_kp,2,dimffnl_kp,ntypat) 333 ! nonlocal form factors at k_prime 334 335 real(dp), pointer :: kinpw_k(:) => null() 336 ! kinpw_k(npw_fft_k) 337 ! (modified) kinetic energy for each plane wave at k 338 ! CAVEAT: In band mode, this array is NOT EQUIVALENT to kinpw(npw_k) 339 340 real(dp), pointer :: kinpw_kp(:) => null() 341 ! kinpw_kp(npw_fft_kp) 342 ! (modified) kinetic energy for each plane wave at k^prime 343 344 real(dp), pointer :: kpg_k(:,:) => null() 345 ! kpg_k(3,npw_fft_k) 346 ! k+G vector coordinates at k 347 348 real(dp), pointer :: kpg_kp(:,:) => null() 349 ! kpg_kp(3,npw_fft_kp) 350 ! k^prime+G vector coordinates at k^prime 351 352 real(dp), ABI_CONTIGUOUS pointer :: phkpxred(:,:) => null() 353 ! phkpxred(2,natom) 354 ! phase factors exp(2 pi k^prime.xred) at k^prime 355 356 real(dp), pointer :: ph3d_k(:,:,:) => null() 357 ! ph3d_k(2,npw_fft_k,matblk) 358 ! 3-dim structure factors, for each atom and plane wave at k 359 360 real(dp), pointer :: ph3d_kp(:,:,:) => null() 361 ! ph3d_kp(2,npw_fft_kp,matblk) 362 ! 3-dim structure factors, for each atom and plane wave at k^prime 363 364 real(dp), pointer :: vectornd(:,:,:,:,:) => null() 365 ! vectornd(n4,n5,n6,nvloc,3) 366 ! vector potential of nuclear magnetic dipoles 367 ! in real space, on the augmented fft grid 368 369 real(dp), pointer :: vlocal(:,:,:,:) => null() 370 ! vlocal(n4,n5,n6,nvloc) 371 ! local potential in real space, on the augmented fft grid 372 373 real(dp), pointer :: vxctaulocal(:,:,:,:,:) => null() 374 ! vxctaulocal(n4,n5,n6,nvloc,4) 375 ! derivative of XC energy density with respect to kinetic energy density, 376 ! in real space, on the augmented fft grid 377 378 real(dp), ABI_CONTIGUOUS pointer :: xred(:,:) => null() 379 ! xred(3,natom) 380 ! reduced coordinates of atoms (dimensionless) 381 382 ! ===== Structured datatype pointers 383 384 type(fock_common_type), pointer :: fockcommon => null() 385 ! common quantities needed to calculate Fock exact exchange 386 387 type(fock_BZ_type), pointer :: fockbz => null() 388 ! total brillouin zone quantities needed to calculate Fock exact exchange 389 390 type(fock_ACE_type), pointer :: fockACE_k => null() 391 ! ACE quantities needed to calculate Fock exact exchange in the ACE context 392 393 contains 394 procedure :: free => destroy_hamiltonian 395 ! Free the memory in the GS Hamiltonian 396 397 procedure :: load_spin => load_spin_hamiltonian 398 ! Setup of the spin-dependent part of the GS Hamiltonian 399 400 procedure :: load_k => load_k_hamiltonian 401 ! Setup of the k-dependent part of the GS Hamiltonian 402 403 procedure :: load_kprime => load_kprime_hamiltonian 404 ! Setup of the k^prime-dependent part of the GS Hamiltonian 405 406 !procedure :: copy => copy_hamiltonian 407 408 end type gs_hamiltonian_type 409 410 public :: init_hamiltonian ! Initialize the GS Hamiltonian 411 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 [fockcommon <type(fock_common_type)>]= common quantities to calculate Fock exact exchange [fockbz <type(fock_BZ_type)>]= quantities to calculate Fock exact exchange in the total BZ 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 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
663 subroutine init_hamiltonian(ham,Psps,pawtab,nspinor,nsppol,nspden,natom,typat,& 664 & xred,nfft,mgfft,ngfft,rprimd,nloalg,& 665 & ph1d,usecprj,comm_atom,mpi_atmtab,mpi_spintab,paw_ij,& ! optional 666 & electronpositron,fock,nucdipmom,use_gpu_cuda) ! optional 667 668 !Arguments ------------------------------------ 669 !scalars 670 integer,intent(in) :: nfft,natom,nspinor,nsppol,nspden,mgfft 671 integer,optional,intent(in) :: comm_atom,usecprj,use_gpu_cuda 672 class(gs_hamiltonian_type),intent(inout),target :: ham 673 type(electronpositron_type),optional,pointer :: electronpositron 674 type(fock_type),optional,pointer :: fock 675 type(pseudopotential_type),intent(in) :: psps 676 !arrays 677 integer,intent(in) :: ngfft(18),nloalg(3),typat(natom) 678 integer,optional,intent(in) :: mpi_atmtab(:),mpi_spintab(2) 679 real(dp),intent(in) :: rprimd(3,3) 680 real(dp),intent(in),target :: xred(3,natom) 681 real(dp),optional,intent(in) :: nucdipmom(3,natom),ph1d(2,3*(2*mgfft+1)*natom) 682 type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw) 683 type(paw_ij_type),optional,intent(in) :: paw_ij(:) 684 685 !Local variables------------------------------- 686 !scalars 687 integer :: my_comm_atom,my_nsppol,itypat,iat,ilmn,indx,isp,cplex_dij,jsp 688 real(dp) :: ucvol 689 !arrays 690 integer :: my_spintab(2) 691 real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3) 692 real(dp),allocatable,target :: ekb_tmp(:,:,:,:) 693 694 ! ************************************************************************* 695 696 DBG_ENTER("COLL") 697 698 !@gs_hamiltonian_type 699 700 !Manage optional parameters 701 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 702 my_spintab=0;my_spintab(1:nsppol)=1;if (present(mpi_spintab)) my_spintab(1:2)=mpi_spintab(1:2) 703 my_nsppol=count(my_spintab==1) 704 705 call metric(gmet,gprimd,-1,rmet,rprimd,ucvol) 706 707 ABI_CHECK(mgfft==MAXVAL(ngfft(1:3)),"Wrong mgfft") 708 709 !Allocate the arrays of the Hamiltonian whose dimensions do not depend on k 710 ABI_MALLOC(ham%atindx,(natom)) 711 ABI_MALLOC(ham%atindx1,(natom)) 712 ABI_MALLOC(ham%typat,(natom)) 713 ABI_MALLOC(ham%indlmn,(6,psps%lmnmax,psps%ntypat)) 714 ABI_MALLOC(ham%nattyp,(psps%ntypat)) 715 ABI_MALLOC(ham%nucdipmom,(3,natom)) 716 ABI_MALLOC(ham%ph1d,(2,3*(2*mgfft+1)*natom)) 717 ABI_MALLOC(ham%pspso,(psps%ntypat)) 718 719 !Initialize most of the Hamiltonian 720 indx=1 721 do itypat=1,psps%ntypat 722 ham%nattyp(itypat)=0 723 do iat=1,natom 724 if (typat(iat)==itypat) then 725 ham%atindx (iat )=indx 726 ham%atindx1(indx)=iat 727 indx=indx+1 728 ham%nattyp(itypat)=ham%nattyp(itypat)+1 729 end if 730 end do 731 end do 732 733 ham%gmet(:,:) =gmet(:,:) 734 ham%gprimd(:,:)=gprimd(:,:) 735 ham%indlmn(:,:,:)=psps%indlmn(:,:,:) 736 ham%lmnmax =psps%lmnmax 737 ham%mgfft =mgfft 738 ham%mpsang =psps%mpsang 739 ham%mpssoang =psps%mpssoang 740 ham%natom =natom 741 ham%nfft =nfft 742 ham%ngfft(:) =ngfft(:) 743 ham%nloalg(:) =nloalg(:) 744 ham%matblk=min(NLO_MINCAT,maxval(ham%nattyp)); if (nloalg(2)>0) ham%matblk=natom 745 ham%nsppol =nsppol 746 ham%nspinor =nspinor 747 ham%ntypat =psps%ntypat 748 ham%typat =typat(1:natom) 749 ham%nvloc=1; if(nspden==4)ham%nvloc=4 750 ham%n4 =ngfft(4) 751 ham%n5 =ngfft(5) 752 ham%n6 =ngfft(6) 753 ham%usepaw =psps%usepaw 754 ham%ucvol =ucvol 755 ham%useylm =psps%useylm 756 ham%use_gpu_cuda=0 ; if(PRESENT(use_gpu_cuda)) ham%use_gpu_cuda=use_gpu_cuda 757 758 ham%pspso(:) =psps%pspso(1:psps%ntypat) 759 if (psps%usepaw==1) then 760 do itypat=1,psps%ntypat 761 ham%pspso(itypat)=1+pawtab(itypat)%usespnorb 762 end do 763 end if 764 765 if (present(nucdipmom)) then 766 ham%nucdipmom(:,:) = nucdipmom(:,:) 767 else 768 ham%nucdipmom(:,:) = zero 769 end if 770 771 ham%xred => xred 772 773 if (present(fock)) then 774 if (associated(fock)) then 775 ham%fockcommon => fock%fock_common 776 if (fock%fock_common%use_ACE==0) ham%fockbz=>fock%fock_BZ 777 end if 778 end if 779 780 if (present(ph1d)) then 781 ham%ph1d(:,:) = ph1d(:,:) 782 else ! Recalculate structure factor phases 783 call getph(ham%atindx,natom,ngfft(1),ngfft(2),ngfft(3),ham%ph1d,xred) 784 end if 785 786 if (ham%usepaw==1) then 787 ham%usecprj=0;if (present(usecprj)) ham%usecprj=usecprj 788 ABI_MALLOC(ham%dimcprj,(natom)) 789 !Be carefull cprj are ordered by atom type (used in non-local operator) 790 call pawcprj_getdim(ham%dimcprj,natom,ham%nattyp,ham%ntypat,ham%typat,pawtab,'O') 791 else 792 ham%usecprj=0 793 ABI_MALLOC(ham%dimcprj,(0)) 794 end if 795 796 ! =========================== 797 ! ==== Non-local factors ==== 798 ! =========================== 799 800 801 if (ham%usepaw==0) then ! Norm-conserving: use constant Kleimann-Bylander energies. 802 ham%dimekb1=psps%dimekb 803 ham%dimekb2=psps%ntypat 804 ham%dimekbq=1 805 ABI_MALLOC(ham%ekb_spin,(psps%dimekb,psps%ntypat,nspinor**2,1,1)) 806 ham%ekb => ham%ekb_spin(:,:,:,:,1) 807 ABI_MALLOC(ham%sij,(0,0)) 808 ham%ekb(:,:,1,1)=psps%ekb(:,:) 809 if (nspinor==2) then 810 ham%ekb(:,:,2,1)=psps%ekb(:,:) 811 ham%ekb(:,:,3:4,1)=zero 812 end if 813 if (PRESENT(electronpositron)) then 814 if (electronpositron_calctype(electronpositron)==1) ham%ekb(:,:,:,:)=-ham%ekb(:,:,:,:) 815 end if 816 817 ! Update enl on GPU (will do it later for PAW) 818 #if defined HAVE_GPU_CUDA 819 if (ham%use_gpu_cuda==1) then 820 call gpu_update_ham_data(ham%ekb(:,:,:,1),size(ham%ekb),ham%sij,size(ham%sij),ham%gprimd,size(ham%gprimd)) 821 end if 822 #endif 823 824 else ! PAW: store overlap coefficients (spin non dependent) and Dij coefficients (spin dependent) 825 cplex_dij=1 826 if (present(paw_ij)) then 827 if (size(paw_ij)>0) cplex_dij=paw_ij(1)%cplex_dij 828 end if 829 if ((nspinor==2).or.any(abs(ham%nucdipmom)>tol8)) cplex_dij=2 830 ham%dimekb1=psps%dimekb*cplex_dij 831 ham%dimekb2=natom 832 ham%dimekbq=1 833 if (present(paw_ij)) then 834 if (size(paw_ij)>0) ham%dimekbq=paw_ij(1)%qphase 835 end if 836 ABI_MALLOC(ham%sij,(ham%dimekb1,psps%ntypat)) 837 do itypat=1,psps%ntypat 838 if (cplex_dij==1) then 839 ham%sij(1:pawtab(itypat)%lmn2_size,itypat)=pawtab(itypat)%sij(:) 840 else 841 do ilmn=1,pawtab(itypat)%lmn2_size 842 ham%sij(2*ilmn-1,itypat)=pawtab(itypat)%sij(ilmn) 843 ham%sij(2*ilmn ,itypat)=zero 844 end do 845 end if 846 if (cplex_dij*pawtab(itypat)%lmn2_size<ham%dimekb1) then 847 ham%sij(cplex_dij*pawtab(itypat)%lmn2_size+1:ham%dimekb1,itypat)=zero 848 end if 849 end do 850 !We preload here PAW non-local factors in order to avoid a communication over atoms 851 ! inside the loop over spins. 852 ABI_MALLOC(ham%ekb_spin,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq,my_nsppol)) 853 ham%ekb_spin=zero 854 if (present(paw_ij)) then 855 if (my_nsppol<ham%nsppol) then 856 ABI_MALLOC(ekb_tmp,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq)) 857 end if 858 jsp=0 859 do isp=1,ham%nsppol 860 if (my_spintab(isp)==1) then 861 jsp=jsp+1 ; ham%ekb => ham%ekb_spin(:,:,:,:,jsp) 862 else 863 ham%ekb => ekb_tmp 864 end if 865 if (present(mpi_atmtab)) then 866 call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom,mpi_atmtab=mpi_atmtab) 867 else 868 call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom) 869 end if 870 end do 871 if (my_nsppol<ham%nsppol) then 872 ABI_FREE(ekb_tmp) 873 end if 874 end if 875 nullify(ham%ekb) 876 end if 877 878 DBG_EXIT("COLL") 879 880 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
1512 subroutine init_rf_hamiltonian(cplex,gs_Ham,ipert,rf_Ham,& 1513 & comm_atom,mpi_atmtab,mpi_spintab,paw_ij1,has_e1kbsc) ! optional arguments 1514 1515 !Arguments ------------------------------------ 1516 !scalars 1517 integer,intent(in) :: cplex,ipert 1518 integer,intent(in),optional :: comm_atom 1519 logical,intent(in),optional :: has_e1kbsc 1520 type(gs_hamiltonian_type),intent(in) :: gs_Ham 1521 type(rf_hamiltonian_type),intent(inout),target :: rf_Ham 1522 !arrays 1523 integer,optional,intent(in) :: mpi_atmtab(:),mpi_spintab(2) 1524 type(paw_ij_type),optional,intent(in) :: paw_ij1(:) 1525 1526 !Local variables------------------------------- 1527 !scalars 1528 integer :: cplex_dij1,isp,jsp,my_comm_atom,my_nsppol 1529 logical :: has_e1kbsc_ 1530 !arrays 1531 integer :: my_spintab(2) 1532 real(dp),allocatable,target :: e1kb_tmp(:,:,:,:) 1533 1534 ! ************************************************************************* 1535 1536 DBG_ENTER("COLL") 1537 1538 !@rf_hamiltonian_type 1539 1540 !Manage optional parameters 1541 has_e1kbsc_=.false.;if (present(has_e1kbsc)) has_e1kbsc_=has_e1kbsc 1542 my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom 1543 my_spintab=0;my_spintab(1:gs_Ham%nsppol)=1;if(present(mpi_spintab)) my_spintab=mpi_spintab 1544 my_nsppol=count(my_spintab==1) 1545 1546 rf_Ham%cplex =cplex 1547 rf_Ham%n4 =gs_Ham%n4 1548 rf_Ham%n5 =gs_Ham%n5 1549 rf_Ham%n6 =gs_Ham%n6 1550 rf_Ham%nvloc =gs_Ham%nvloc 1551 rf_Ham%nsppol =gs_Ham%nsppol 1552 rf_Ham%nspinor =gs_Ham%nspinor 1553 1554 rf_Ham%dime1kb1=0 1555 rf_Ham%dime1kb2=gs_Ham%dimekb2 1556 if (gs_Ham%usepaw==1.and.ipert/=gs_Ham%natom+1.and.ipert/=gs_Ham%natom+10) then 1557 cplex_dij1=1;if ((gs_Ham%nspinor==2).or.any(abs(gs_Ham%nucdipmom)>tol8)) cplex_dij1=2 1558 rf_Ham%dime1kb1=cplex_dij1*(gs_Ham%lmnmax*(gs_Ham%lmnmax+1))/2 1559 end if 1560 1561 ! Allocate the arrays of the 1st-order Hamiltonian 1562 ! We preload here 1st-order non-local factors in order to avoid 1563 ! a communication over atoms inside the loop over spins. 1564 if (gs_Ham%usepaw==1.and.rf_Ham%dime1kb1>0) then 1565 if ((ipert>=1.and.ipert<=gs_Ham%natom).or.ipert==gs_Ham%natom+2.or.& 1566 ipert==gs_Ham%natom+3.or.ipert==gs_Ham%natom+4.or.ipert==gs_Ham%natom+11) then 1567 1568 ABI_MALLOC(rf_Ham%e1kbfr_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol)) 1569 rf_Ham%e1kbfr_spin=zero 1570 if (has_e1kbsc_) then 1571 ABI_MALLOC(rf_Ham%e1kbsc_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol)) 1572 rf_Ham%e1kbsc_spin=zero 1573 end if 1574 1575 if (present(paw_ij1)) then 1576 1577 if (my_nsppol<rf_Ham%nsppol) then 1578 ABI_MALLOC(e1kb_tmp,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex)) 1579 end if 1580 1581 ! === Frozen term 1582 jsp=0 1583 do isp=1,rf_Ham%nsppol 1584 if (my_spintab(isp)==1) then 1585 jsp=jsp+1 ; rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsp) 1586 else 1587 rf_Ham%e1kbfr => e1kb_tmp 1588 end if 1589 if (present(mpi_atmtab)) then 1590 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr,mpi_atmtab=mpi_atmtab) 1591 else 1592 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr) 1593 end if 1594 end do 1595 1596 ! === Self-consistent term 1597 if (has_e1kbsc_) then 1598 jsp=0 1599 do isp=1,rf_Ham%nsppol 1600 if (my_spintab(isp)==1) then 1601 jsp=jsp+1 ; rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsp) 1602 else 1603 rf_Ham%e1kbsc => e1kb_tmp 1604 end if 1605 if (present(mpi_atmtab)) then 1606 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc,mpi_atmtab=mpi_atmtab) 1607 else 1608 call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc) 1609 end if 1610 end do 1611 end if 1612 1613 if (my_nsppol<rf_Ham%nsppol) then 1614 ABI_FREE(e1kb_tmp) 1615 end if 1616 1617 end if 1618 end if 1619 end if 1620 1621 if (.not.allocated(rf_Ham%e1kbfr_spin)) then 1622 ABI_MALLOC(rf_Ham%e1kbfr_spin,(0,0,0,0,0)) 1623 end if 1624 if (.not.allocated(rf_Ham%e1kbsc_spin)) then 1625 ABI_MALLOC(rf_Ham%e1kbsc_spin,(0,0,0,0,0)) 1626 end if 1627 nullify(rf_Ham%e1kbfr) 1628 nullify(rf_Ham%e1kbsc) 1629 1630 DBG_EXIT("COLL") 1631 1632 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
916 subroutine load_k_hamiltonian(ham,ffnl_k,fockACE_k,gbound_k,istwf_k,kinpw_k,& 917 kg_k,kpg_k,kpt_k,npw_k,npw_fft_k,ph3d_k,& 918 compute_gbound,compute_ph3d) 919 920 !Arguments ------------------------------------ 921 !scalars 922 integer,intent(in),optional :: npw_k,npw_fft_k,istwf_k 923 logical,intent(in),optional :: compute_gbound,compute_ph3d 924 class(gs_hamiltonian_type),intent(inout),target :: ham 925 !arrays 926 integer,intent(in),optional,target :: gbound_k(:,:),kg_k(:,:) 927 real(dp),intent(in),optional :: kpt_k(3) 928 real(dp),intent(in),optional,target :: ffnl_k(:,:,:,:),kinpw_k(:),kpg_k(:,:),ph3d_k(:,:,:) 929 type(fock_ACE_type),intent(in),optional,target :: fockACE_k 930 931 !Local variables------------------------------- 932 !scalars 933 integer :: iat,iatom 934 logical :: compute_gbound_ 935 real(dp) :: arg 936 !character(len=500) :: msg 937 938 ! ************************************************************************* 939 940 DBG_ENTER("COLL") 941 942 !@gs_hamiltonian_type 943 944 !k-dependent scalars 945 if (present(kpt_k)) then 946 ham%kpt_k(:) = kpt_k(:) 947 ham%kpt_kp(:) = kpt_k(:) 948 end if 949 if (present(istwf_k)) then 950 ham%istwf_k = istwf_k 951 ham%istwf_kp = istwf_k 952 end if 953 if (present(npw_k)) then 954 ham%npw_k = npw_k 955 ham%npw_kp = npw_k 956 end if 957 if (present(npw_fft_k)) then 958 ham%npw_fft_k = npw_fft_k 959 ham%npw_fft_kp = npw_fft_k 960 else if (present(npw_k)) then 961 ham%npw_fft_k = npw_k 962 ham%npw_fft_kp = npw_k 963 end if 964 965 !Pointers to k-dependent quantitites 966 if (present(kinpw_k)) then 967 ham%kinpw_k => kinpw_k 968 ham%kinpw_kp => kinpw_k 969 end if 970 if (present(kg_k)) then 971 ham%kg_k => kg_k 972 ham%kg_kp => kg_k 973 end if 974 if (present(kpg_k)) then 975 ham%kpg_k => kpg_k 976 ham%kpg_kp => kpg_k 977 end if 978 if (present(ffnl_k)) then 979 ham%ffnl_k => ffnl_k 980 ham%ffnl_kp => ffnl_k 981 end if 982 if (present(ph3d_k)) then 983 ham%ph3d_k => ph3d_k 984 ham%ph3d_kp => ph3d_k 985 end if 986 if (present(fockACE_k)) then 987 ham%fockACE_k => fockACE_k 988 end if 989 !Compute exp(i.k.R) for each atom 990 if (present(kpt_k)) then 991 if (associated(Ham%phkpxred).and.(.not.associated(Ham%phkpxred,Ham%phkxred))) then 992 ABI_FREE(Ham%phkpxred) 993 end if 994 ABI_SFREE(ham%phkxred) 995 ABI_MALLOC(ham%phkxred,(2,ham%natom)) 996 do iat=1,ham%natom 997 iatom=ham%atindx(iat) 998 arg=two_pi*DOT_PRODUCT(kpt_k,ham%xred(:,iat)) 999 ham%phkxred(1,iatom)=DCOS(arg) 1000 ham%phkxred(2,iatom)=DSIN(arg) 1001 end do 1002 ham%phkpxred => ham%phkxred 1003 end if 1004 1005 !Compute or copy G sphere boundary at k+g 1006 compute_gbound_=.false.;if (present(compute_gbound)) compute_gbound_=compute_gbound 1007 if (present(gbound_k)) compute_gbound_=.true. 1008 if (compute_gbound_) then 1009 if (associated(Ham%gbound_kp,Ham%gbound_k)) then 1010 nullify(Ham%gbound_kp) 1011 else if (associated(Ham%gbound_kp)) then 1012 ABI_FREE(Ham%gbound_kp) 1013 end if 1014 ABI_SFREE(ham%gbound_k) 1015 end if 1016 if (.not.allocated(ham%gbound_k)) then 1017 ABI_MALLOC(ham%gbound_k,(2*ham%mgfft+8,2)) 1018 ham%gbound_k(:,:)=0 1019 ham%gbound_kp => ham%gbound_k 1020 end if 1021 if (compute_gbound_) then 1022 if (present(gbound_k)) then 1023 ham%gbound_k(:,:)=gbound_k(:,:) 1024 else 1025 if (.not.associated(ham%kg_k)) then 1026 ABI_BUG('Something is missing for gbound_k computation!') 1027 end if 1028 !write(std_out,*)"About to call sphereboundary" 1029 !write(std_out,*)"size(kg_k), npw_k, mgfft",size(ham%kg_k, dim=2), ham%npw_k, ham%mgfft 1030 call sphereboundary(ham%gbound_k,ham%istwf_k,ham%kg_k,ham%mgfft,ham%npw_k) 1031 end if 1032 ham%gbound_kp => ham%gbound_k 1033 end if 1034 1035 !Compute 3D structure factors for each atom at k+g 1036 if (present(compute_ph3d).and.present(ph3d_k)) then 1037 if (compute_ph3d.and.ham%nloalg(2)>0) then 1038 if ((.not.allocated(ham%phkxred)).or.(.not.associated(ham%kg_k)).or.& 1039 (.not.associated(ham%ph3d_k))) then 1040 ABI_BUG('Something is missing for ph3d_k computation!') 1041 end if 1042 call ph1d3d(1,ham%natom,ham%kg_k,ham%matblk,ham%natom,ham%npw_k,ham%ngfft(1),& 1043 ham%ngfft(2),ham%ngfft(3),ham%phkxred,ham%ph1d,ham%ph3d_k) 1044 end if 1045 end if 1046 1047 DBG_EXIT("COLL") 1048 1049 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
1727 subroutine load_k_rf_hamiltonian(rf_Ham,dkinpw_k,ddkinpw_k,npw_k) 1728 1729 !Arguments ------------------------------------ 1730 !scalars 1731 integer,intent(in),optional :: npw_k 1732 class(rf_hamiltonian_type),intent(inout),target :: rf_Ham 1733 !arrays 1734 real(dp),intent(in),optional,target :: dkinpw_k(:),ddkinpw_k(:) 1735 1736 !Local variables------------------------------- 1737 1738 ! ************************************************************************* 1739 1740 DBG_ENTER("COLL") 1741 1742 !@gs_hamiltonian_type 1743 1744 !k-dependent scalars 1745 if (present(npw_k)) then 1746 rf_Ham%npw_k = npw_k 1747 rf_Ham%npw_kp = npw_k 1748 end if 1749 1750 !Pointers to k-dependent quantitites 1751 if (present(dkinpw_k)) then 1752 rf_Ham%dkinpw_k => dkinpw_k 1753 rf_Ham%dkinpw_kp => dkinpw_k 1754 end if 1755 if (present(ddkinpw_k)) then 1756 rf_Ham%ddkinpw_k => ddkinpw_k 1757 rf_Ham%ddkinpw_kp => ddkinpw_k 1758 end if 1759 1760 DBG_EXIT("COLL") 1761 1762 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
1084 subroutine load_kprime_hamiltonian(ham,ffnl_kp,gbound_kp,istwf_kp,kinpw_kp,& 1085 kg_kp,kpg_kp,kpt_kp,npw_kp,npw_fft_kp,& 1086 ph3d_kp,compute_gbound,compute_ph3d) 1087 1088 !Arguments ------------------------------------ 1089 !scalars 1090 integer,intent(in),optional :: npw_kp,npw_fft_kp,istwf_kp 1091 logical,intent(in),optional :: compute_gbound,compute_ph3d 1092 class(gs_hamiltonian_type),intent(inout),target :: ham 1093 !arrays 1094 integer,intent(in),optional,target :: gbound_kp(:,:),kg_kp(:,:) 1095 real(dp),intent(in),optional :: kpt_kp(3) 1096 real(dp),intent(in),optional,target :: ffnl_kp(:,:,:,:),kinpw_kp(:),kpg_kp(:,:),ph3d_kp(:,:,:) 1097 1098 !Local variables------------------------------- 1099 !scalars 1100 integer :: iat,iatom 1101 logical :: compute_gbound_ 1102 real(dp) :: arg 1103 !character(len=500) :: msg 1104 1105 ! ************************************************************************* 1106 1107 DBG_ENTER("COLL") 1108 1109 !@gs_hamiltonian_type 1110 1111 !k-dependent scalars 1112 if (present(kpt_kp)) ham%kpt_kp(:)= kpt_kp(:) 1113 if (present(istwf_kp)) ham%istwf_kp = istwf_kp 1114 if (present(npw_kp)) ham%npw_kp = npw_kp 1115 if (present(npw_fft_kp)) then 1116 ham%npw_fft_kp = npw_fft_kp 1117 else if (present(npw_kp)) then 1118 ham%npw_fft_kp = npw_kp 1119 end if 1120 1121 !Pointers to k-dependent quantitites 1122 if (present(kinpw_kp)) ham%kinpw_kp => kinpw_kp 1123 if (present(kg_kp)) ham%kg_kp => kg_kp 1124 if (present(kpg_kp)) ham%kpg_kp => kpg_kp 1125 if (present(ffnl_kp)) ham%ffnl_kp => ffnl_kp 1126 if (present(ph3d_kp)) ham%ph3d_kp => ph3d_kp 1127 1128 !Compute exp(i.k^prime.R) for each atom 1129 if (present(kpt_kp)) then 1130 if (associated(ham%phkpxred,ham%phkxred)) then 1131 nullify(ham%phkpxred) 1132 else if (associated(ham%phkpxred)) then 1133 ABI_FREE(ham%phkpxred) 1134 end if 1135 ABI_MALLOC(ham%phkpxred,(2,ham%natom)) 1136 do iat=1,ham%natom 1137 iatom=ham%atindx(iat) 1138 arg=two_pi*DOT_PRODUCT(kpt_kp,ham%xred(:,iat)) 1139 ham%phkpxred(1,iatom)=DCOS(arg) 1140 ham%phkpxred(2,iatom)=DSIN(arg) 1141 end do 1142 end if 1143 1144 !Compute or copy G sphere boundary at k^prime+g 1145 compute_gbound_=.false. 1146 if (present(kpt_kp).and.present(compute_gbound)) compute_gbound_=compute_gbound 1147 if (present(gbound_kp)) compute_gbound_=.true. 1148 if (compute_gbound_) then 1149 if (associated(ham%gbound_kp,ham%gbound_k)) then 1150 nullify(ham%gbound_kp) 1151 else if (associated(ham%gbound_kp)) then 1152 ABI_FREE(ham%gbound_kp) 1153 end if 1154 if (present(gbound_kp)) then 1155 ham%gbound_kp(:,:)=gbound_kp(:,:) 1156 else 1157 if (.not.associated(ham%kg_kp)) then 1158 ABI_BUG('Something is missing for gbound_kp computation!') 1159 end if 1160 ABI_MALLOC(ham%gbound_kp,(2*ham%mgfft+8,2)) 1161 call sphereboundary(ham%gbound_kp,ham%istwf_kp,ham%kg_kp,ham%mgfft,ham%npw_kp) 1162 end if 1163 end if 1164 1165 !Compute 3D structure factors for each atom at k^prime+g 1166 if (present(compute_ph3d).and.present(ph3d_kp)) then 1167 if (compute_ph3d.and.ham%nloalg(2)>0) then 1168 if ((.not.associated(ham%phkpxred)).or.(.not.associated(ham%kg_kp)).or.& 1169 & (.not.associated(ham%ph3d_kp))) then 1170 ABI_BUG('Something is missing for ph3d_kp computation!') 1171 end if 1172 call ph1d3d(1,ham%natom,ham%kg_kp,ham%matblk,ham%natom,ham%npw_kp,ham%ngfft(1),& 1173 & ham%ngfft(2),ham%ngfft(3),ham%phkpxred,ham%ph1d,ham%ph3d_kp) 1174 end if 1175 end if 1176 1177 DBG_EXIT("COLL") 1178 1179 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
1382 subroutine load_spin_hamiltonian(Ham,isppol,vectornd,vlocal,vxctaulocal,with_nonlocal) 1383 1384 !Arguments ------------------------------------ 1385 !scalars 1386 integer,intent(in) :: isppol 1387 logical,optional,intent(in) :: with_nonlocal 1388 class(gs_hamiltonian_type),intent(inout),target :: Ham 1389 !arrays 1390 real(dp),optional,intent(in),target :: vectornd(:,:,:,:,:) 1391 real(dp),optional,intent(in),target :: vlocal(:,:,:,:),vxctaulocal(:,:,:,:,:) 1392 1393 !Local variables------------------------------- 1394 !scalars 1395 integer :: jsppol 1396 1397 ! ************************************************************************* 1398 1399 DBG_ENTER("COLL") 1400 1401 !@gs_hamiltonian_type 1402 if (present(vlocal)) then 1403 ABI_CHECK(size(vlocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc, "Wrong vlocal") 1404 Ham%vlocal => vlocal 1405 end if 1406 if (present(vxctaulocal)) then 1407 ABI_CHECK(size(vxctaulocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*4, "Wrong vxctaulocal") 1408 Ham%vxctaulocal => vxctaulocal 1409 end if 1410 if (present(vectornd)) then 1411 ABI_CHECK(size(vectornd)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*3, "Wrong vectornd") 1412 Ham%vectornd => vectornd 1413 end if 1414 1415 ! Retrieve non-local factors for this spin component 1416 if (present(with_nonlocal)) then 1417 if (with_nonlocal) then 1418 jsppol=min(isppol,size(Ham%ekb_spin,5)) 1419 if (jsppol>0) Ham%ekb => Ham%ekb_spin(:,:,:,:,jsppol) 1420 end if 1421 end if 1422 1423 ! Update enl and sij on GPU 1424 #if defined HAVE_GPU_CUDA 1425 if (Ham%use_gpu_cuda==1) then 1426 call gpu_update_ham_data(Ham%ekb(:,:,:,1),size(Ham%ekb),Ham%sij,size(Ham%sij),Ham%gprimd,size(Ham%gprimd)) 1427 end if 1428 #endif 1429 1430 DBG_EXIT("COLL") 1431 1432 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 [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
1657 subroutine load_spin_rf_hamiltonian(rf_Ham,isppol,vectornd,vlocal1,with_nonlocal) 1658 1659 !Arguments ------------------------------------ 1660 !scalars 1661 integer,intent(in) :: isppol 1662 logical,optional,intent(in) :: with_nonlocal 1663 class(rf_hamiltonian_type),intent(inout),target :: rf_Ham 1664 !arrays 1665 real(dp),optional,target,intent(in) :: vlocal1(:,:,:,:) 1666 real(dp),optional,target,intent(in) :: vectornd(:,:,:,:) 1667 1668 !Local variables------------------------------- 1669 !scalars 1670 integer :: jsppol 1671 1672 ! ************************************************************************* 1673 1674 DBG_ENTER("COLL") 1675 1676 !@rf_hamiltonian_type 1677 1678 if (present(vlocal1)) then 1679 ABI_CHECK(size(vlocal1)==rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc,"Wrong vlocal1") 1680 rf_Ham%vlocal1 => vlocal1 1681 end if 1682 1683 if (present(vectornd)) then 1684 ABI_CHECK(size(vectornd)==rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc,"Wrong vectornd") 1685 rf_Ham%vectornd => vectornd 1686 end if 1687 1688 ! Retrieve non-local factors for this spin component 1689 if (present(with_nonlocal)) then 1690 if (with_nonlocal) then 1691 if (size(rf_Ham%e1kbfr_spin)>0) then 1692 jsppol=min(isppol,size(rf_Ham%e1kbfr_spin,5)) 1693 if (jsppol>0) rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsppol) 1694 end if 1695 if (size(rf_Ham%e1kbsc_spin)>0) then 1696 jsppol=min(isppol,size(rf_Ham%e1kbsc_spin,5)) 1697 if (jsppol>0) rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsppol) 1698 end if 1699 end if 1700 end if 1701 1702 DBG_EXIT("COLL") 1703 1704 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
1859 subroutine pawdij2e1kb(paw_ij1,isppol,comm_atom,mpi_atmtab,e1kbfr,e1kbsc) 1860 1861 !Arguments ------------------------------------ 1862 !scalars 1863 integer,intent(in) :: isppol,comm_atom 1864 !arrays 1865 integer,intent(in),optional,target :: mpi_atmtab(:) 1866 real(dp),optional,intent(out) :: e1kbfr(:,:,:,:),e1kbsc(:,:,:,:) 1867 type(paw_ij_type),intent(in) :: paw_ij1(:) 1868 1869 !Local variables------------------------------- 1870 !scalars 1871 integer :: dimdij1,dime1kb1,dime1kb3,dime1kb4,iatom,iatom_tot,ierr,isp,ispden 1872 integer :: my_natom,natom,qphase 1873 logical :: my_atmtab_allocated,paral_atom 1874 !arrays 1875 integer,pointer :: my_atmtab(:) 1876 1877 ! ************************************************************************* 1878 1879 DBG_ENTER("COLL") 1880 1881 if ((.not.present(e1kbfr)).and.(.not.present(e1kbsc))) return 1882 if (present(e1kbfr)) then 1883 e1kbfr=zero ; natom=size(e1kbfr,2) 1884 end if 1885 if (present(e1kbsc)) then 1886 e1kbsc=zero ; natom=size(e1kbsc,2) 1887 end if 1888 1889 !Set up parallelism over atoms 1890 my_natom=size(paw_ij1) ; paral_atom=(xmpi_comm_size(comm_atom)>1) 1891 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 1892 call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 1893 1894 !Retrieve 1st-order PAW Dij coefficients for this spin component (frozen) 1895 if (my_natom>0.and.present(e1kbfr)) then 1896 if (allocated(paw_ij1(1)%dijfr)) then 1897 dime1kb1=size(e1kbfr,1) ; dime1kb3=size(e1kbfr,3) ; dime1kb4=size(e1kbfr,4) 1898 ABI_CHECK(paw_ij1(1)%qphase==dime1kb4,'BUG in pawdij2e1kb (1)!') 1899 do ispden=1,dime1kb3 1900 isp=isppol;if (dime1kb3==4) isp=ispden 1901 do iatom=1,my_natom 1902 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 1903 qphase=paw_ij1(iatom)%qphase 1904 dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size 1905 ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!') 1906 e1kbfr(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dijfr(1:dimdij1,isp) 1907 if (qphase==2) e1kbfr(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp) 1908 end do 1909 end do 1910 end if 1911 end if 1912 1913 !Retrieve 1st-order PAW Dij coefficients for this spin component (self-consistent) 1914 if (my_natom>0.and.present(e1kbsc)) then 1915 if (allocated(paw_ij1(1)%dijfr).and.allocated(paw_ij1(1)%dij)) then 1916 dime1kb1=size(e1kbsc,1) ; dime1kb3=size(e1kbsc,3) ; dime1kb4=size(e1kbsc,4) 1917 ABI_CHECK(paw_ij1(1)%qphase==dime1kb4,'BUG in pawdij2e1kb (1)!') 1918 do ispden=1,dime1kb3 1919 isp=isppol;if (dime1kb3==4) isp=ispden 1920 do iatom=1,my_natom 1921 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 1922 qphase=paw_ij1(iatom)%qphase 1923 dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size 1924 ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!') 1925 e1kbsc(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dij (1:dimdij1,isp) & 1926 & -paw_ij1(iatom)%dijfr(1:dimdij1,isp) 1927 if (qphase==2) e1kbsc(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dij (dimdij1+1:2*dimdij1,isp) & 1928 & -paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp) 1929 end do 1930 end do 1931 end if 1932 end if 1933 1934 ! Communication in case of distribution over atomic sites 1935 if (paral_atom) then 1936 if (present(e1kbfr)) call xmpi_sum(e1kbfr,comm_atom,ierr) 1937 if (present(e1kbsc)) call xmpi_sum(e1kbsc,comm_atom,ierr) 1938 end if 1939 1940 !Destroy atom table used for parallelism 1941 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 1942 1943 DBG_EXIT("COLL") 1944 1945 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
1781 subroutine pawdij2ekb(ekb,paw_ij,isppol,comm_atom,mpi_atmtab) 1782 1783 !Arguments ------------------------------------ 1784 !scalars 1785 integer,intent(in) :: isppol,comm_atom 1786 !arrays 1787 integer,intent(in),optional,target :: mpi_atmtab(:) 1788 real(dp),intent(out) :: ekb(:,:,:,:) 1789 type(paw_ij_type),intent(in) :: paw_ij(:) 1790 1791 !Local variables------------------------------- 1792 !scalars 1793 integer :: dimdij,dimekb1,dimekb3,dimekb4,iatom,iatom_tot,ierr,ii,isp,ispden,my_natom,natom,qphase 1794 logical :: my_atmtab_allocated,paral_atom 1795 !arrays 1796 integer,pointer :: my_atmtab(:) 1797 1798 ! ************************************************************************* 1799 1800 DBG_ENTER("COLL") 1801 1802 ekb=zero 1803 1804 !Set up parallelism over atoms 1805 natom=size(ekb,2); my_natom=size(paw_ij) 1806 paral_atom=(xmpi_comm_size(comm_atom)>1) 1807 nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab 1808 call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom) 1809 1810 !Retrieve PAW Dij coefficients for this spin component 1811 if (my_natom>0) then 1812 if (allocated(paw_ij(1)%dij)) then 1813 dimekb1=size(ekb,1) ; dimekb3=size(ekb,3) ; dimekb4=size(ekb,4) 1814 qphase=paw_ij(1)%qphase 1815 ABI_CHECK(qphase<=dimekb4,'paw_ij%qphase>dimekb4!') 1816 do ii=1,qphase 1817 do ispden=1,dimekb3 1818 isp=isppol; if (dimekb3==4) isp=ispden 1819 do iatom=1,my_natom 1820 iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom) 1821 dimdij=paw_ij(iatom)%cplex_dij*paw_ij(iatom)%lmn2_size 1822 ABI_CHECK(dimdij<=dimekb1,'Size of paw_ij%dij>dimekb1!') 1823 ekb(1:dimdij,iatom_tot,ispden,ii)=paw_ij(iatom)%dij(1+(ii-1)*dimdij:ii*dimdij,isp) 1824 end do 1825 end do 1826 end do 1827 end if 1828 end if 1829 1830 !Communication in case of distribution over atomic sites 1831 if (paral_atom) then 1832 call xmpi_sum(ekb,comm_atom,ierr) 1833 end if 1834 1835 !Destroy atom table used for parallelism 1836 call free_my_atmtab(my_atmtab,my_atmtab_allocated) 1837 1838 DBG_EXIT("COLL") 1839 1840 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
427 type,public :: rf_hamiltonian_type 428 429 ! ===== Integer scalars 430 431 integer :: cplex 432 ! if 1, real space 1-order functions on FFT grid are REAL; if 2, COMPLEX 433 434 integer :: dime1kb1 435 ! First dimension of E1kb, derivative of Ekb with respect to a perturbation 436 437 integer :: dime1kb2 438 ! Second dimension of E1kb, derivative of Ekb with respect to a perturbation 439 ! NCPP: dime1kb2=ntypat, PAW: dime1kb2=natom 440 441 integer :: npw_k 442 ! number of plane waves at k 443 444 integer :: npw_kp 445 ! number of plane waves at k^prime 446 447 integer:: nspinor 448 ! Number of spinorial components 449 450 integer :: nsppol 451 ! Total number of spin components (1=non-polarized, 2=polarized) 452 453 integer :: nvloc 454 ! Number of components of vloc 455 ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4 456 457 integer :: n4,n5,n6 458 ! same as ngfft(4:6) 459 460 ! ===== Real arrays 461 462 real(dp), allocatable :: e1kbfr_spin(:,:,:,:,:) 463 ! e1kbfr_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol) 464 ! Contains the values of e1kbfr array for all spins treated by current process 465 ! See e1kbfr description ; e1kbfr is pointer to e1kbfr_spin(:,:,:,:,isppol) 466 467 real(dp), allocatable :: e1kbsc_spin(:,:,:,:,:) 468 ! e1kbsc_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol) 469 ! Contains the values of e1kbsc array for all spins treated by current process 470 ! See e1kbsc description ; e1kbsc is pointer to e1kbsc_spin(:,:,:,:,isppol) 471 472 ! ===== Real pointers 473 474 real(dp), pointer :: dkinpw_k(:) => null() 475 ! dkinpw_k(npw_k) 476 ! 1st derivative of the (modified) kinetic energy for each plane wave at k 477 478 real(dp), pointer :: dkinpw_kp(:) => null() 479 ! dkinpw_kp(npw_kp) 480 ! 1st derivative of the (modified) kinetic energy for each plane wave at k^prime 481 482 real(dp), pointer :: ddkinpw_k(:) => null() 483 ! ddkinpw_k(npw_k) 484 ! 2nd derivative of the (modified) kinetic energy for each plane wave at k 485 486 real(dp), pointer :: ddkinpw_kp(:) => null() 487 ! ddkinpw_kp(npw_kp) 488 ! 2nd derivative of the (modified) kinetic energy for each plane wave at k^prime 489 490 real(dp), pointer :: e1kbfr(:,:,:,:) => null() 491 ! Frozen part of 1st derivative of ekb for the considered perturbation 492 ! (part not depending on VHxc^(1)) 493 ! e1kbfr(dime1kb1,dime1kb2,nspinor**2,cplex) 494 ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol) 495 496 real(dp), ABI_CONTIGUOUS pointer :: e1kbsc(:,:,:,:) => null() 497 ! Self-consistent 1st derivative of ekb for the considered perturbation 498 ! (part depending only on self-consistent VHxc^(1)) 499 ! e1kbsc(dime1kb1,dime1kb2,nspinor**2,cplex) 500 ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol) 501 502 real(dp), pointer :: vectornd(:,:,:,:) => null() 503 ! vectornd(n4,n5,n6,nvloc) 504 ! vector potential of nuclear magnetic dipoles 505 ! in real space, on the augmented fft grid, in direction idir 506 ! (the ddk pert direction) 507 508 real(dp), pointer :: vlocal1(:,:,:,:) => null() 509 ! vlocal1(cplex*n4,n5,n6,nvloc) 510 ! 1st-order local potential in real space, on the augmented fft grid 511 512 contains 513 procedure :: free => destroy_rf_hamiltonian 514 ! Free the memory in the RF Hamiltonian 515 516 procedure :: load_spin => load_spin_rf_hamiltonian 517 ! Setup of the spin-dependent part of the RF Hamiltonian. 518 519 procedure :: load_k => load_k_rf_hamiltonian 520 ! Setup of the k-dependent part of the RF Hamiltonian 521 522 end type rf_hamiltonian_type 523 524 public :: init_rf_hamiltonian ! Initialize the RF Hamiltonian