TABLE OF CONTENTS
ABINIT/prep_bandfft_tabs [ Functions ]
NAME
prep_bandfft_tabs
FUNCTION
This routine transpose various tabs needed in bandfft parallelization
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (FB,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 . For the initials of contributors, see ~abinit/doc/developers/contributors.txt.
INPUTS
ikpt=index of the k-point gs_hamk <type(gs_hamiltonian_type)>=all data for the Hamiltonian at k mkmem=number of k points which can fit in memory; set to 0 if use disk
SIDE EFFECTS
bandfft_kpt tabs (defined in m_bandfft_kpt module)
PARENTS
energy,fock2ACE,forstrnps,vtorho
CHILDREN
bandfft_kpt_init2,bandfft_kpt_set_ikpt,mkkpg,timab,xmpi_allgatherv
NOTES
SOURCE
34 #if defined HAVE_CONFIG_H 35 #include "config.h" 36 #endif 37 38 #include "abi_common.h" 39 40 subroutine prep_bandfft_tabs(gs_hamk,ikpt,mkmem,mpi_enreg) 41 42 use defs_basis 43 use defs_abitypes 44 use m_errors 45 use m_profiling_abi 46 use m_xmpi 47 use m_bandfft_kpt, only : bandfft_kpt,bandfft_kpt_set_ikpt,bandfft_kpt_init2 48 use m_hamiltonian, only : gs_hamiltonian_type 49 50 !This section has been created automatically by the script Abilint (TD). 51 !Do not modify the following lines by hand. 52 #undef ABI_FUNC 53 #define ABI_FUNC 'prep_bandfft_tabs' 54 use interfaces_18_timing 55 use interfaces_66_nonlocal 56 !End of the abilint section 57 58 implicit none 59 60 !Arguments ------------------------------- 61 integer,intent(in) :: ikpt,mkmem 62 type(gs_hamiltonian_type),intent(inout) :: gs_hamk 63 type(MPI_type),intent(inout) :: mpi_enreg 64 65 !Local variables------------------------------- 66 integer :: dimffnl,ierr,ikpt_this_proc,ipw,lmnmax,matblk,ndatarecv,nkpg,npw_k,ntypat,spaceComm 67 logical :: tabs_allocated 68 real(dp) :: tsec(2) 69 character(len=500) :: message 70 integer, allocatable :: recvcounts(:),rdispls(:) 71 integer, allocatable :: recvcountsloc(:),rdisplsloc(:) 72 real(dp),allocatable :: ffnl_gather(:,:,:,:),ffnl_little(:,:,:,:),ffnl_little_gather(:,:,:,:) 73 real(dp),allocatable :: kinpw_gather(:),kpg_k_gather(:,:) 74 real(dp),allocatable :: ph3d_gather(:,:,:),ph3d_little(:,:,:),ph3d_little_gather(:,:,:) 75 76 ! ********************************************************************* 77 78 call timab(575,1,tsec) 79 80 ikpt_this_proc=mpi_enreg%my_kpttab(ikpt) 81 tabs_allocated=((bandfft_kpt(ikpt_this_proc)%flag1_is_allocated==1).and.& 82 & (ikpt_this_proc <= mkmem).and.(ikpt_this_proc/=0)) 83 84 if (.not.tabs_allocated) then 85 message = ' the bandfft tabs are not allocated !' 86 MSG_BUG(message) 87 end if 88 89 ntypat=gs_hamk%ntypat 90 lmnmax=gs_hamk%lmnmax 91 matblk=gs_hamk%matblk 92 npw_k=gs_hamk%npw_k 93 dimffnl=size(gs_hamk%ffnl_k,2) 94 nkpg=size(gs_hamk%kpg_k,2) 95 96 ABI_ALLOCATE(rdispls,(mpi_enreg%nproc_band)) 97 ABI_ALLOCATE(recvcounts,(mpi_enreg%nproc_band)) 98 ABI_ALLOCATE(rdisplsloc,(mpi_enreg%nproc_band)) 99 ABI_ALLOCATE(recvcountsloc,(mpi_enreg%nproc_band)) 100 101 spaceComm =mpi_enreg%comm_band 102 ndatarecv =bandfft_kpt(ikpt_this_proc)%ndatarecv 103 rdispls(:) =bandfft_kpt(ikpt_this_proc)%rdispls(:) 104 recvcounts(:)=bandfft_kpt(ikpt_this_proc)%recvcounts(:) 105 106 !---- Process FFNL 107 if (associated(gs_hamk%ffnl_k)) then 108 ABI_ALLOCATE(ffnl_gather,(ndatarecv,dimffnl,lmnmax,ntypat)) 109 ABI_ALLOCATE(ffnl_little,(dimffnl,lmnmax,ntypat,npw_k)) 110 ABI_ALLOCATE(ffnl_little_gather,(dimffnl,lmnmax,ntypat,ndatarecv)) 111 do ipw=1,npw_k 112 ffnl_little(:,:,:,ipw)=gs_hamk%ffnl_k(ipw,:,:,:) 113 end do 114 recvcountsloc(:)=recvcounts(:)*dimffnl*lmnmax*ntypat 115 rdisplsloc(:)=rdispls(:)*dimffnl*lmnmax*ntypat 116 call xmpi_allgatherv(ffnl_little,npw_k*dimffnl*lmnmax*ntypat,ffnl_little_gather,& 117 & recvcountsloc(:),rdisplsloc,spaceComm,ierr) 118 do ipw=1,ndatarecv 119 ffnl_gather(ipw,:,:,:)=ffnl_little_gather(:,:,:,ipw) 120 end do 121 ABI_DEALLOCATE(ffnl_little) 122 ABI_DEALLOCATE(ffnl_little_gather) 123 else 124 ABI_ALLOCATE(ffnl_gather,(0,0,0,0)) 125 end if 126 127 !---- Process PH3D 128 if (associated(gs_hamk%ph3d_k)) then 129 ABI_ALLOCATE(ph3d_gather,(2,ndatarecv,matblk)) 130 ABI_ALLOCATE(ph3d_little,(2,matblk,npw_k)) 131 ABI_ALLOCATE(ph3d_little_gather,(2,matblk,ndatarecv)) 132 recvcountsloc(:)=recvcounts(:)*2*matblk 133 rdisplsloc(:)=rdispls(:)*2*matblk 134 do ipw=1,npw_k 135 ph3d_little(:,:,ipw)=gs_hamk%ph3d_k(:,ipw,:) 136 end do 137 call xmpi_allgatherv(ph3d_little,npw_k*2*matblk,ph3d_little_gather,& 138 & recvcountsloc(:),rdisplsloc,spaceComm,ierr) 139 do ipw=1,ndatarecv 140 ph3d_gather(:,ipw,:)=ph3d_little_gather(:,:,ipw) 141 end do 142 ABI_DEALLOCATE(ph3d_little_gather) 143 ABI_DEALLOCATE(ph3d_little) 144 else 145 ABI_ALLOCATE(ph3d_gather,(0,0,0)) 146 end if 147 148 !---- Process KPG_K 149 if (associated(gs_hamk%kpg_k)) then 150 ABI_ALLOCATE(kpg_k_gather,(ndatarecv,nkpg)) 151 if (nkpg>0) then 152 call mkkpg(bandfft_kpt(ikpt_this_proc)%kg_k_gather,kpg_k_gather,gs_hamk%kpt_k,nkpg,ndatarecv) 153 end if 154 else 155 ABI_ALLOCATE(kpg_k_gather,(0,0)) 156 end if 157 158 !---- Process KINPW 159 if (associated(gs_hamk%kinpw_k)) then 160 ABI_ALLOCATE(kinpw_gather,(ndatarecv)) 161 recvcountsloc(:)=recvcounts(:) 162 rdisplsloc(:)=rdispls(:) 163 call xmpi_allgatherv(gs_hamk%kinpw_k,npw_k,kinpw_gather,recvcountsloc(:),rdisplsloc,spaceComm,ierr) 164 else 165 ABI_ALLOCATE(kinpw_gather,(0)) 166 end if 167 168 ABI_DEALLOCATE(recvcounts) 169 ABI_DEALLOCATE(rdispls) 170 ABI_DEALLOCATE(recvcountsloc) 171 ABI_DEALLOCATE(rdisplsloc) 172 173 call bandfft_kpt_init2(bandfft_kpt,dimffnl,ffnl_gather,ikpt_this_proc,kinpw_gather,kpg_k_gather,& 174 & lmnmax,matblk,mkmem,ndatarecv,nkpg,ntypat,ph3d_gather) 175 176 ABI_DEALLOCATE(ffnl_gather) 177 ABI_DEALLOCATE(ph3d_gather) 178 ABI_DEALLOCATE(kinpw_gather) 179 ABI_DEALLOCATE(kpg_k_gather) 180 181 !---- Store current kpt index 182 call bandfft_kpt_set_ikpt(ikpt,mpi_enreg) 183 184 call timab(575,2,tsec) 185 186 end subroutine prep_bandfft_tabs