TABLE OF CONTENTS
- ABINIT/m_rttddft_exponential
- m_rttddft_exponential/paral_kgb_transpose
- m_rttddft_exponential/rttddft_exp_taylor
ABINIT/m_rttddft_exponential [ Modules ]
NAME
m_rttddft_exponential
FUNCTION
Contains subroutines to compute the exponential part of the propagator using various approximations
COPYRIGHT
Copyright (C) 2021-2022 ABINIT group (FB) 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
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 module m_rttddft_exponential 24 25 use defs_basis 26 use defs_abitypes, only: MPI_type 27 28 use m_bandfft_kpt, only: bandfft_kpt, & 29 & bandfft_kpt_get_ikpt 30 use m_dtset, only: dataset_type 31 use m_getghc, only: multithreaded_getghc 32 use m_hamiltonian, only: gs_hamiltonian_type 33 use m_invovl, only: apply_invovl 34 use m_pawcprj, only: pawcprj_type, pawcprj_alloc, pawcprj_free 35 use m_prep_kgb, only: prep_getghc, prep_index_wavef_bandpp 36 use m_profiling_abi, only: abimem_record 37 use m_rttddft_properties, only: rttddft_calc_eig, rttddft_calc_enl 38 use m_xmpi, only: xmpi_alltoallv 39 40 implicit none 41 42 private
m_rttddft_exponential/paral_kgb_transpose [ Functions ]
[ Top ] [ m_rttddft_exponential ] [ Functions ]
NAME
paral_kgb_transpose
FUNCTION
if option = 1: Forward transpose Transpose cg_1 in linalg ((npw/npband),nband) distribution into cg_2 in fft (npw,bandpp) distribution if option = -1: Backward transpose Transpose back cg_2 in fft (npw,bandpp) distribution into cg_1 linakg (npw/npband),nband) distribution
INPUTS
if option = 1: cg_1 <real((npw/nband)*nspinor*nband)> if option = -1: cg_2 <real(npw*nspinor*bandpp)> nband_t <integer> = number of bands after forward transpose (bandpp) npw_t <integer> = number of pw after forward transpose (npw_k) dtset <type(dataset_type)> = all input variables for this dataset mpi_enreg <MPI_type> = MPI-parallelisation information nspinor <integer> = "dimension of spinors" option <integer> = option for forward or backward transpose index_wavef_band <integer> = order of the bands after transpose
OUTPUT
if option = 1 : cg_2 <real(npw*nspinor*bandpp)> nband_t <integer> = number of bands after forward transpose (bandpp) npw_t <integer> = number of pw after forward transpose (npw_k) if option = -1: cg_1 <real((npw/nband)*nspinor*nband)>
SIDE EFFECTS
if option = 1 : cg_2 has been allocated index_wavef_band has been allocated if option = -1: cg_2 has been deallocated index_wavef_band has been allocated
SOURCE
268 subroutine paral_kgb_transpose(cg_1,cg_2,mpi_enreg,nband_t,npw_t,nspinor,option,index_wavef_band) 269 270 implicit none 271 272 !Arguments ------------------------------------ 273 !scalars 274 integer, intent(inout) :: nband_t 275 integer, intent(inout) :: npw_t 276 integer, intent(in) :: nspinor 277 integer, intent(in) :: option 278 type(MPI_type), intent(inout) :: mpi_enreg 279 !arrays 280 integer, allocatable, intent(inout) :: index_wavef_band(:) 281 real(dp), intent(inout) :: cg_1(:,:) 282 real(dp), allocatable, intent(inout) :: cg_2(:,:) 283 284 !Local variables------------------------------- 285 !scalars 286 integer :: bandpp 287 integer :: ierr 288 integer :: ikpt_this_proc 289 !arrays 290 integer :: recvcountsloc(mpi_enreg%nproc_band) 291 integer :: rdisplsloc(mpi_enreg%nproc_band) 292 integer :: sendcountsloc(mpi_enreg%nproc_band) 293 integer :: sdisplsloc(mpi_enreg%nproc_band) 294 real(dp), allocatable :: cg_work(:,:) 295 296 ! *********************************************************************** 297 298 !Init useful MPI variables 299 bandpp = mpi_enreg%bandpp 300 ikpt_this_proc = bandfft_kpt_get_ikpt() 301 recvcountsloc = bandfft_kpt(ikpt_this_proc)%recvcounts*2*nspinor*bandpp 302 rdisplsloc = bandfft_kpt(ikpt_this_proc)%rdispls*2*nspinor*bandpp 303 sendcountsloc = bandfft_kpt(ikpt_this_proc)%sendcounts*2*nspinor 304 sdisplsloc = bandfft_kpt(ikpt_this_proc)%sdispls*2*nspinor 305 306 !Forward transpose: cg_1 -> cg_2 307 if (option == 1) then 308 nband_t = bandpp 309 npw_t = bandfft_kpt(ikpt_this_proc)%ndatarecv 310 ABI_MALLOC(cg_2, (2, npw_t*nspinor*nband_t)) 311 ABI_MALLOC(cg_work, (2, npw_t*nspinor*nband_t)) 312 !Transpose input cg_1 into cg_work 313 call xmpi_alltoallv(cg_1,sendcountsloc,sdisplsloc,cg_work, & 314 & recvcountsloc,rdisplsloc,mpi_enreg%comm_band,ierr) 315 !properly sort array according to bandd after alltoall 316 call prep_index_wavef_bandpp(mpi_enreg%nproc_band,mpi_enreg%bandpp, & 317 & nspinor,bandfft_kpt(ikpt_this_proc)%ndatarecv , & 318 & bandfft_kpt(ikpt_this_proc)%recvcounts, & 319 & bandfft_kpt(ikpt_this_proc)%rdispls, index_wavef_band) 320 cg_2(:,:) = cg_work(:,index_wavef_band) 321 ABI_FREE(cg_work) 322 end if 323 324 !Transpose back: cg_2 -> cg_1 325 if (option == -1) then 326 ABI_MALLOC(cg_work, (2, npw_t*nspinor*nband_t)) 327 cg_work(:,index_wavef_band) = cg_2(:,:) 328 !Transpose cg_work to input cg_1 329 call xmpi_alltoallv(cg_work,recvcountsloc,rdisplsloc,cg_1, & 330 & sendcountsloc,sdisplsloc,mpi_enreg%comm_band,ierr) 331 !Free memory 332 ABI_FREE(cg_work) 333 ABI_FREE(cg_2) 334 ABI_FREE(index_wavef_band) 335 end if 336 337 end subroutine paral_kgb_transpose
m_rttddft_exponential/rttddft_exp_taylor [ Functions ]
[ Top ] [ m_rttddft_exponential ] [ Functions ]
NAME
rttddft_propagator_exp
FUNCTION
Applies the propagator exp(-i*dt*H(*S^-1)) on the cg coeffs approximating the exponential using a Taylor exapansion
INPUTS
cg <real(npw*nspinor*nband)> = the wavefunction coefficients dtset <type(dataset_type)>=all input variables for this dataset ham_k <type(gs_hamiltonian_type)> = the Hamiltonian H mpi_enreg <MPI_type> = MPI-parallelisation information nband_k <integer> = number of bands npw_k <integer> = number of plane waves nspinor <integer> = dimension of spinors
OUTPUT
cg <real(npw*nspinor*nband)> = the new cg after application of the exponential propagator enl <real(bandpp)> = non local contribution to the energy in the NC case - optional eig <real(bandpp)> = eigenvalues - optional
SOURCE
77 subroutine rttddft_exp_taylor(cg,dtset,ham_k,mpi_enreg,nband_k,npw_k,nspinor,enl,eig) 78 79 implicit none 80 81 !Arguments ------------------------------------ 82 !scalars 83 integer, intent(in) :: nband_k 84 integer, intent(in) :: npw_k 85 integer, intent(in) :: nspinor 86 type(dataset_type), intent(in) :: dtset 87 type(gs_hamiltonian_type), intent(inout) :: ham_k 88 type(MPI_type), intent(inout) :: mpi_enreg 89 !arrays 90 real(dp), target, intent(inout) :: cg(2,npw_k*nband_k*nspinor) 91 real(dp), intent(out), optional :: enl(:) 92 real(dp), pointer, intent(inout), optional :: eig(:) 93 94 !Local variables------------------------------- 95 !scalars 96 integer :: cpopt 97 integer :: iorder 98 integer :: nband_t 99 integer :: npw_t 100 integer :: sij_opt 101 integer :: tim_getghc = 5 102 integer :: nfact 103 logical :: l_paw, l_eig, l_enl 104 real(dp) :: dt 105 !arrays 106 integer, allocatable :: index_wavef_band(:) 107 type(pawcprj_type), allocatable :: cwaveprj(:,:) 108 real(dp), pointer :: cg_t(:,:) => null() 109 real(dp), target, allocatable :: cg_trans(:,:) 110 real(dp), allocatable :: ghc(:,:) 111 real(dp), allocatable :: gvnlxc_dummy(:,:) 112 real(dp), allocatable :: gsm1hc(:,:) 113 real(dp), allocatable :: gsc(:,:) 114 real(dp), allocatable :: tmp(:,:) 115 116 ! *********************************************************************** 117 118 !Check what properties we need to compute 119 l_eig = .false.; l_enl = .false. 120 if (present(eig)) then 121 if (associated(eig)) l_eig = .true. 122 end if 123 if (present(enl)) l_enl = (size(enl)>0) 124 125 !Transpose if paral_kgb 126 if (dtset%paral_kgb == 1 .and. mpi_enreg%nproc_band>1) then 127 call paral_kgb_transpose(cg,cg_trans,mpi_enreg,nband_t,npw_t,nspinor,1,index_wavef_band) 128 cg_t => cg_trans 129 else 130 npw_t = npw_k 131 nband_t = nband_k 132 cg_t => cg 133 end if 134 135 l_paw = (ham_k%usepaw == 1) 136 if(l_paw) then 137 ABI_MALLOC(cwaveprj, (ham_k%natom,nspinor*nband_t)) 138 call pawcprj_alloc(cwaveprj,0,ham_k%dimcprj) 139 else 140 ABI_MALLOC(cwaveprj,(0,0)) 141 end if 142 cpopt = 0 143 144 dt = dtset%dtele !electronic timestep 145 146 ABI_MALLOC(ghc, (2, npw_t*nspinor*nband_t)) 147 ABI_MALLOC(gsc, (2, npw_t*nspinor*nband_t)) 148 ABI_MALLOC(gsm1hc,(2, npw_t*nspinor*nband_t)) 149 ABI_MALLOC(gvnlxc_dummy, (0, 0)) 150 151 !*** Taylor expansion *** 152 ABI_MALLOC(tmp,(2, npw_t*nspinor*nband_t)) 153 tmp(:,:) = cg_t(:,:) 154 nfact = 1 155 156 do iorder = 1, dtset%td_exp_order 157 158 nfact = nfact*iorder 159 160 !** Apply Hamiltonian 161 sij_opt = 0 162 if (iorder == 1 .and. l_eig .and. l_paw) sij_opt = 1 163 if (dtset%paral_kgb /= 1) then 164 call multithreaded_getghc(cpopt,tmp,cwaveprj,ghc,gsc,ham_k,gvnlxc_dummy,1.0_dp, & 165 & mpi_enreg,nband_t,dtset%prtvol,sij_opt,tim_getghc,0) 166 else 167 call prep_getghc(tmp,ham_k,gvnlxc_dummy,ghc,gsc,1.0_dp,nband_t,mpi_enreg, & 168 & dtset%prtvol,sij_opt,cpopt,cwaveprj,already_transposed=.true.) 169 end if 170 171 !** Also apply S^-1 in PAW case 172 if (l_paw) then 173 call apply_invovl(ham_k,ghc,gsm1hc,cwaveprj,npw_t,nband_t,mpi_enreg,nspinor,dtset%diago_apply_block_sliced) 174 tmp(1,:) = dt*gsm1hc(2,:)/real(nfact,dp) 175 tmp(2,:) = -dt*gsm1hc(1,:)/real(nfact,dp) 176 else 177 tmp(1,:) = dt*ghc(2,:)/real(nfact,dp) 178 tmp(2,:) = -dt*ghc(1,:)/real(nfact,dp) 179 end if 180 181 !** Compute properties if requested 182 if (iorder == 1) then 183 !eigenvalues 184 if (l_eig) then 185 if (l_paw) then 186 call rttddft_calc_eig(cg_t,eig,ghc,ham_k%istwf_k,nband_t,npw_t,nspinor, & 187 & mpi_enreg%me_g0,mpi_enreg%comm_spinorfft,gsc=gsc) 188 else 189 call rttddft_calc_eig(cg_t,eig,ghc,ham_k%istwf_k,nband_t,npw_t,nspinor, & 190 & mpi_enreg%me_g0,mpi_enreg%comm_spinorfft) 191 end if 192 end if 193 194 !non local PSP energy contribution in NC case 195 if (l_enl) then 196 call rttddft_calc_enl(cg_t,enl,ham_k,nband_t,npw_t,nspinor,mpi_enreg) 197 end if 198 end if 199 200 !** Update cg 201 cg_t(:,:) = cg_t(:,:) + tmp(:,:) 202 203 end do 204 205 ABI_FREE(tmp) 206 ABI_FREE(ghc) 207 ABI_FREE(gsc) 208 ABI_FREE(gsm1hc) 209 ABI_FREE(gvnlxc_dummy) 210 211 if(l_paw) call pawcprj_free(cwaveprj) 212 ABI_FREE(cwaveprj) 213 214 !Transpose back if paral_kgb 215 if (dtset%paral_kgb == 1 .and. mpi_enreg%nproc_band > 1) then 216 call paral_kgb_transpose(cg,cg_trans,mpi_enreg,nband_t,npw_t,nspinor,-1,index_wavef_band) 217 end if 218 219 nullify(cg_t) 220 221 end subroutine rttddft_exp_taylor