Branch data Line data Source code
1 : : !!****m* ABINIT/m_rttddft_exponential
2 : : !! NAME
3 : : !! m_rttddft_exponential
4 : : !!
5 : : !! FUNCTION
6 : : !! Contains subroutines to compute the exponential
7 : : !! part of the propagator using various approximations
8 : : !!
9 : : !! COPYRIGHT
10 : : !! Copyright (C) 2021-2026 ABINIT group (FB)
11 : : !! This file is distributed under the terms of the
12 : : !! GNU General Public License, see ~abinit/COPYING
13 : : !! or http://www.gnu.org/copyleft/gpl.txt .
14 : : !!
15 : : !! SOURCE
16 : :
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
43 : : !!***
44 : :
45 : : public :: rttddft_exp_taylor
46 : : !!***
47 : :
48 : : contains
49 : : !!***
50 : :
51 : : !!****f* m_rttddft_exponential/rttddft_exp_taylor
52 : : !!
53 : : !! NAME
54 : : !! rttddft_propagator_exp
55 : : !!
56 : : !! FUNCTION
57 : : !! Applies the propagator exp(-i*dt*H(*S^-1)) on the cg coeffs
58 : : !! approximating the exponential using a Taylor exapansion
59 : : !!
60 : : !! INPUTS
61 : : !! cg <real(npw*nspinor*nband)> = the wavefunction coefficients
62 : : !! dtset <type(dataset_type)>=all input variables for this dataset
63 : : !! ham_k <type(gs_hamiltonian_type)> = the Hamiltonian H
64 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
65 : : !! nband_k <integer> = number of bands
66 : : !! npw_k <integer> = number of plane waves
67 : : !! nspinor <integer> = dimension of spinors
68 : : !!
69 : : !! OUTPUT
70 : : !! cg <real(npw*nspinor*nband)> = the new cg after application
71 : : !! of the exponential propagator
72 : : !! enl <real(bandpp)> = non local contribution to the energy in the
73 : : !! NC case - optional
74 : : !! eig <real(bandpp)> = eigenvalues - optional
75 : : !!
76 : : !! SOURCE
77 [ + + + - : 27380 : subroutine rttddft_exp_taylor(cg,dtset,ham_k,mpi_enreg,nband_k,npw_k,nspinor,enl,eig)
- + ]
78 : :
79 : : !Arguments ------------------------------------
80 : : !scalars
81 : : integer, intent(in) :: nband_k
82 : : integer, intent(in) :: npw_k
83 : : integer, intent(in) :: nspinor
84 : : type(dataset_type), intent(in) :: dtset
85 : : type(gs_hamiltonian_type), intent(inout) :: ham_k
86 : : type(MPI_type), intent(inout) :: mpi_enreg
87 : : !arrays
88 : : real(dp), target, intent(inout) :: cg(2,npw_k*nband_k*nspinor)
89 : : real(dp), intent(out), optional :: enl(:)
90 : : real(dp), pointer, intent(inout), optional :: eig(:)
91 : :
92 : : !Local variables-------------------------------
93 : : !scalars
94 : : integer :: cpopt
95 : : integer :: iorder
96 : : integer :: nband_t
97 : : integer :: npw_t
98 : : integer :: sij_opt
99 : : integer :: tim_getghc = 5
100 : : logical :: l_paw, l_eig, l_enl
101 : : real(dp) :: dt
102 : : !arrays
103 : 27380 : integer, allocatable :: index_wavef_band(:)
104 : 27380 : type(pawcprj_type), allocatable :: cwaveprj(:,:)
105 : : real(dp), pointer :: cg_t(:,:) => null()
106 : 27380 : real(dp), target, allocatable :: cg_trans(:,:)
107 : 27380 : real(dp), allocatable :: ghc(:,:)
108 : 27380 : real(dp), allocatable :: gvnlxc_dummy(:,:)
109 : 27380 : real(dp), allocatable :: gsm1hc(:,:)
110 : 27380 : real(dp), allocatable :: gsc(:,:)
111 : 27380 : real(dp), allocatable :: tmp(:,:)
112 : :
113 : : ! ***********************************************************************
114 : :
115 : : !Check what properties we need to compute
116 : 27380 : l_eig = .false.; l_enl = .false.
117 [ + + ]: 27380 : if (present(eig)) then
118 [ + - ]: 27260 : if (associated(eig)) l_eig = .true.
119 : : end if
120 [ + + ]: 27380 : if (present(enl)) l_enl = (size(enl)>0)
121 : :
122 : : !Transpose if paral_kgb
123 [ + + + - ]: 27380 : if (dtset%paral_kgb == 1 .and. mpi_enreg%nproc_band>1) then
124 : 6560 : call paral_kgb_transpose(cg,cg_trans,mpi_enreg,nband_t,npw_t,nspinor,1,index_wavef_band)
125 : 6560 : cg_t => cg_trans
126 : : else
127 : 20820 : npw_t = npw_k
128 : 20820 : nband_t = nband_k
129 : 20820 : cg_t => cg
130 : : end if
131 : :
132 : 27380 : l_paw = (ham_k%usepaw == 1)
133 [ + + ]: 27380 : if(l_paw) then
134 [ - + + - : 724400 : ABI_MALLOC(cwaveprj, (ham_k%natom,nspinor*nband_t))
- + + - +
- + - - +
- + - + +
+ + + +
- ]
135 : 27020 : call pawcprj_alloc(cwaveprj,0,ham_k%dimcprj)
136 : : else
137 [ - + - + ]: 360 : ABI_MALLOC(cwaveprj,(0,0))
138 : : end if
139 : 27380 : cpopt = 0
140 : :
141 : 27380 : dt = dtset%dtele !electronic timestep
142 : :
143 [ - + + - : 82140 : ABI_MALLOC(ghc, (2, npw_t*nspinor*nband_t))
+ - + - -
+ - + ]
144 [ - + + - : 54760 : ABI_MALLOC(gsc, (2, npw_t*nspinor*nband_t))
- + - + -
+ ]
145 [ - + + - : 54760 : ABI_MALLOC(gsm1hc,(2, npw_t*nspinor*nband_t))
- + - + ]
146 [ - + ]: 27380 : ABI_MALLOC(gvnlxc_dummy, (0, 0))
147 : :
148 : : !*** Taylor expansion ***
149 [ - + + - : 54760 : ABI_MALLOC(tmp,(2, npw_t*nspinor*nband_t))
- + - + ]
150 [ + + + + ]: 158948300 : tmp(:,:) = cg_t(:,:)
151 : :
152 [ + + ]: 136900 : do iorder = 1, dtset%td_exp_order
153 : :
154 : : !** Apply Hamiltonian
155 : 109520 : sij_opt = 0
156 [ - + - - ]: 109520 : if (iorder == 1 .and. l_eig .and. l_paw) sij_opt = 1
157 [ + + ]: 109520 : if (dtset%paral_kgb /= 1) then
158 : : call multithreaded_getghc(cpopt,tmp,cwaveprj,ghc,gsc,ham_k,gvnlxc_dummy,1.0_dp, &
159 : 83280 : & mpi_enreg,nband_t,dtset%prtvol,sij_opt,tim_getghc,0)
160 : : else
161 : : call prep_getghc(tmp,ham_k,gvnlxc_dummy,ghc,gsc,1.0_dp,nband_t,mpi_enreg, &
162 : 26240 : & dtset%prtvol,sij_opt,cpopt,cwaveprj,already_transposed=.true.)
163 : : end if
164 : :
165 : : !** Also apply S^-1 in PAW case
166 [ + + ]: 109520 : if (l_paw) then
167 : 108080 : call apply_invovl(ham_k,ghc,gsm1hc,cwaveprj,npw_t,nband_t,mpi_enreg,nspinor,dtset%invovl_blksliced,dtset%invovl_precision)
168 [ + + ]: 176904080 : tmp(1,:) = dt*gsm1hc(2,:)/real(iorder,dp)
169 [ + + ]: 176904080 : tmp(2,:) = -dt*gsm1hc(1,:)/real(iorder,dp)
170 : : else
171 [ + + ]: 35100000 : tmp(1,:) = dt*ghc(2,:)/real(iorder,dp)
172 [ + + ]: 35100000 : tmp(2,:) = -dt*ghc(1,:)/real(iorder,dp)
173 : : end if
174 : :
175 : : !** Compute properties if requested
176 [ + + ]: 109520 : if (iorder == 1) then
177 : : !eigenvalues
178 [ - + ]: 27380 : if (l_eig) then
179 [ # # ]: 0 : if (l_paw) then
180 : : call rttddft_calc_eig(cg_t,eig,ghc,ham_k%istwf_k,nband_t,npw_t,nspinor, &
181 [ # # # # ]: 0 : & mpi_enreg%me_g0,mpi_enreg%comm_spinorfft,gsc=gsc)
182 : : else
183 : : call rttddft_calc_eig(cg_t,eig,ghc,ham_k%istwf_k,nband_t,npw_t,nspinor, &
184 [ # # # # ]: 0 : & mpi_enreg%me_g0,mpi_enreg%comm_spinorfft)
185 : : end if
186 : : end if
187 : :
188 : : !non local PSP energy contribution in NC case
189 [ + + ]: 27380 : if (l_enl) then
190 [ - + - - : 320 : call rttddft_calc_enl(cg_t,enl,ham_k,nband_t,npw_t,nspinor,mpi_enreg)
- + - + -
- ]
191 : : end if
192 : : end if
193 : :
194 : : !** Update cg
195 [ + + + + ]: 635820580 : cg_t(:,:) = cg_t(:,:) + tmp(:,:)
196 : :
197 : : end do
198 : :
199 [ - + ]: 27380 : ABI_FREE(tmp)
200 [ - + ]: 27380 : ABI_FREE(ghc)
201 [ - + ]: 27380 : ABI_FREE(gsc)
202 [ - + ]: 27380 : ABI_FREE(gsm1hc)
203 [ - + ]: 27380 : ABI_FREE(gvnlxc_dummy)
204 : :
205 [ + + ]: 27380 : if(l_paw) call pawcprj_free(cwaveprj)
206 [ + - + + : 444980 : ABI_FREE(cwaveprj)
- + - + -
+ ]
207 : :
208 : : !Transpose back if paral_kgb
209 [ + + + - ]: 27380 : if (dtset%paral_kgb == 1 .and. mpi_enreg%nproc_band > 1) then
210 : 6560 : call paral_kgb_transpose(cg,cg_trans,mpi_enreg,nband_t,npw_t,nspinor,-1,index_wavef_band)
211 : : end if
212 : :
213 : 27380 : nullify(cg_t)
214 : :
215 [ - + - + : 27380 : end subroutine rttddft_exp_taylor
- + - - -
- - - -
+ ]
216 : : !!***
217 : :
218 : : !!****f* m_rttddft_exponential/paral_kgb_transpose
219 : : !!
220 : : !! NAME
221 : : !! paral_kgb_transpose
222 : : !!
223 : : !! FUNCTION
224 : : !! if option = 1: Forward transpose
225 : : !! Transpose cg_1 in linalg ((npw/npband),nband) distribution
226 : : !! into cg_2 in fft (npw,bandpp) distribution
227 : : !!
228 : : !! if option = -1: Backward transpose
229 : : !! Transpose back cg_2 in fft (npw,bandpp) distribution
230 : : !! into cg_1 linakg (npw/npband),nband) distribution
231 : : !!
232 : : !! INPUTS
233 : : !! if option = 1:
234 : : !! cg_1 <real((npw/nband)*nspinor*nband)>
235 : : !! if option = -1:
236 : : !! cg_2 <real(npw*nspinor*bandpp)>
237 : : !! nband_t <integer> = number of bands after forward transpose (bandpp)
238 : : !! npw_t <integer> = number of pw after forward transpose (npw_k)
239 : : !! dtset <type(dataset_type)> = all input variables for this dataset
240 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
241 : : !! nspinor <integer> = "dimension of spinors"
242 : : !! option <integer> = option for forward or backward transpose
243 : : !! index_wavef_band <integer> = order of the bands after transpose
244 : : !!
245 : : !! OUTPUT
246 : : !! if option = 1 :
247 : : !! cg_2 <real(npw*nspinor*bandpp)>
248 : : !! nband_t <integer> = number of bands after forward transpose (bandpp)
249 : : !! npw_t <integer> = number of pw after forward transpose (npw_k)
250 : : !! if option = -1:
251 : : !! cg_1 <real((npw/nband)*nspinor*nband)>
252 : :
253 : : !! SIDE EFFECTS
254 : : !! if option = 1 :
255 : : !! cg_2 has been allocated
256 : : !! index_wavef_band has been allocated
257 : : !! if option = -1:
258 : : !! cg_2 has been deallocated
259 : : !! index_wavef_band has been allocated
260 : : !!
261 : : !! SOURCE
262 [ - + ]: 13120 : subroutine paral_kgb_transpose(cg_1,cg_2,mpi_enreg,nband_t,npw_t,nspinor,option,index_wavef_band)
263 : :
264 : : !Arguments ------------------------------------
265 : : !scalars
266 : : integer, intent(inout) :: nband_t
267 : : integer, intent(inout) :: npw_t
268 : : integer, intent(in) :: nspinor
269 : : integer, intent(in) :: option
270 : : type(MPI_type), intent(inout) :: mpi_enreg
271 : : !arrays
272 : : integer, allocatable, intent(inout) :: index_wavef_band(:)
273 : : real(dp), intent(inout) :: cg_1(:,:)
274 : : real(dp), allocatable, intent(inout) :: cg_2(:,:)
275 : :
276 : : !Local variables-------------------------------
277 : : !scalars
278 : : integer :: bandpp
279 : : integer :: ierr
280 : : integer :: ikpt_this_proc
281 : : !arrays
282 : 26240 : integer :: recvcountsloc(mpi_enreg%nproc_band)
283 : 26240 : integer :: rdisplsloc(mpi_enreg%nproc_band)
284 : 26240 : integer :: sendcountsloc(mpi_enreg%nproc_band)
285 : 26240 : integer :: sdisplsloc(mpi_enreg%nproc_band)
286 : 13120 : real(dp), allocatable :: cg_work(:,:)
287 : :
288 : : ! ***********************************************************************
289 : :
290 : : !Init useful MPI variables
291 : 13120 : bandpp = mpi_enreg%bandpp
292 : 13120 : ikpt_this_proc = bandfft_kpt_get_ikpt()
293 [ + + ]: 39360 : recvcountsloc = bandfft_kpt(ikpt_this_proc)%recvcounts*2*nspinor*bandpp
294 [ + + ]: 39360 : rdisplsloc = bandfft_kpt(ikpt_this_proc)%rdispls*2*nspinor*bandpp
295 [ + + ]: 39360 : sendcountsloc = bandfft_kpt(ikpt_this_proc)%sendcounts*2*nspinor
296 [ + + ]: 39360 : sdisplsloc = bandfft_kpt(ikpt_this_proc)%sdispls*2*nspinor
297 : :
298 : : !Forward transpose: cg_1 -> cg_2
299 [ + + ]: 13120 : if (option == 1) then
300 : 6560 : nband_t = bandpp
301 : 6560 : npw_t = bandfft_kpt(ikpt_this_proc)%ndatarecv
302 [ - + + - : 19680 : ABI_MALLOC(cg_2, (2, npw_t*nspinor*nband_t))
+ - + - -
+ - + -
+ ]
303 [ - + + - : 13120 : ABI_MALLOC(cg_work, (2, npw_t*nspinor*nband_t))
- + - + ]
304 : : !Transpose input cg_1 into cg_work
305 : : call xmpi_alltoallv(cg_1,sendcountsloc,sdisplsloc,cg_work, &
306 : 6560 : & recvcountsloc,rdisplsloc,mpi_enreg%comm_band,ierr)
307 : : !properly sort array according to bandd after alltoall
308 : : call prep_index_wavef_bandpp(mpi_enreg%nproc_band,mpi_enreg%bandpp, &
309 : : & nspinor,bandfft_kpt(ikpt_this_proc)%ndatarecv , &
310 : : & bandfft_kpt(ikpt_this_proc)%recvcounts, &
311 [ - + ]: 6560 : & bandfft_kpt(ikpt_this_proc)%rdispls, index_wavef_band)
312 [ + + + + ]: 18759680 : cg_2(:,:) = cg_work(:,index_wavef_band)
313 [ - + ]: 6560 : ABI_FREE(cg_work)
314 : : end if
315 : :
316 : : !Transpose back: cg_2 -> cg_1
317 [ + + ]: 13120 : if (option == -1) then
318 [ - + + - : 19680 : ABI_MALLOC(cg_work, (2, npw_t*nspinor*nband_t))
+ - + - -
+ - + -
+ ]
319 [ + + + + ]: 18759680 : cg_work(:,index_wavef_band) = cg_2(:,:)
320 : : !Transpose cg_work to input cg_1
321 : : call xmpi_alltoallv(cg_work,recvcountsloc,rdisplsloc,cg_1, &
322 : 6560 : & sendcountsloc,sdisplsloc,mpi_enreg%comm_band,ierr)
323 : : !Free memory
324 [ - + ]: 6560 : ABI_FREE(cg_work)
325 [ - + ]: 6560 : ABI_FREE(cg_2)
326 [ - + ]: 6560 : ABI_FREE(index_wavef_band)
327 : : end if
328 : :
329 [ - + ]: 13120 : end subroutine paral_kgb_transpose
330 : : !!***
331 : :
332 : : end module m_rttddft_exponential
333 : : !!***
|
