Branch data Line data Source code
1 : : !!****m* ABINIT/m_fock_getghc
2 : : !! NAME
3 : : !!
4 : : !! FUNCTION
5 : : !!
6 : : !! COPYRIGHT
7 : : !! Copyright (C) 2013-2026 ABINIT group (CMartins, FJ, MT, XG)
8 : : !! This file is distributed under the terms of the
9 : : !! GNU General Public License, see ~abinit/COPYING
10 : : !! or http://www.gnu.org/copyleft/gpl.txt .
11 : : !!
12 : : !! SOURCE
13 : :
14 : : #if defined HAVE_CONFIG_H
15 : : #include "config.h"
16 : : #endif
17 : :
18 : : #include "abi_common.h"
19 : :
20 : : module m_fock_getghc
21 : :
22 : : use, intrinsic :: iso_c_binding, only: c_size_t, c_loc
23 : : use defs_basis
24 : : use m_abicore
25 : : use m_errors
26 : : use m_xmpi
27 : : use m_fock
28 : : use m_pawcprj
29 : : !use m_cgtools
30 : :
31 : : use defs_abitypes, only : mpi_type
32 : : use defs_datatypes, only : pseudopotential_type
33 : : use m_time, only : timab, time_accu
34 : : use m_matrix, only : matr3inv
35 : : use m_cgtools, only : dotprod_g,dotprod_g_batch_half,dotprod_g_batch_full
36 : : use m_kg, only : mkkpg
37 : : use m_fftcore, only : sphereboundary
38 : : use m_fft, only : fftpac, fourwf, fourdp
39 : : use m_fstrings, only : sjoin, itoa
40 : : use m_hamiltonian, only : gs_hamiltonian_type, K_H_KPRIME
41 : : use m_paw_nhat, only : pawmknhat_psipsi_ndat, pawdijhat_ndat
42 : : use m_spacepar, only : hartre
43 : : use m_nonlop, only : nonlop
44 : : use m_bandfft_kpt, only : bandfft_kpt, bandfft_kpt_type, bandfft_kpt_savetabs,bandfft_kpt_restoretabs, &
45 : : prep_bandfft_tabs
46 : : use m_pawtab, only : pawtab_type
47 : : use m_paw_ij, only : paw_ij_type
48 : : use m_mkffnl, only : mkffnl
49 : : use m_mpinfo, only : proc_distrb_cycle
50 : : use m_gputk
51 : : use m_abi_linalg
52 : :
53 : : #if defined(HAVE_GPU)
54 : : use m_gpu_toolbox
55 : : #endif
56 : :
57 : : implicit none
58 : :
59 : : private
60 : : !!***
61 : :
62 : : public :: fock_getghc
63 : : public :: fock2ACE
64 : : public :: fock_ACE_getghc
65 : : !!***
66 : :
67 : : contains
68 : : !!***
69 : :
70 : 0 : subroutine select_ndat_occ_for_gpu(ndat_occ,nband_k,ndat,npw,cplex_fock,nfftf,ngfft,&
71 : 0 : n4,n5,n6,natom,nspinor,lmn2_size,usepaw,cprj,ieigen,need_ghc,optfor,optstr,nfgd_max)
72 : :
73 : : !Arguments ------------------------------------
74 : : ! Scalars
75 : : integer,intent(in) :: nband_k,ndat,npw,cplex_fock,nfftf,n4,n5,n6,ngfft(18)
76 : : integer,intent(in) :: natom,nspinor,lmn2_size,usepaw,ieigen,nfgd_max
77 : : logical,intent(in) :: optfor,optstr,need_ghc
78 : : integer,intent(out) :: ndat_occ
79 : : type(pawcprj_type),intent(in) :: cprj(natom,nspinor*ndat)
80 : :
81 : : !Local variables-------------------------------
82 : : integer :: i,ider,nprojs
83 : : integer,target :: t_fft(3)
84 : : integer(kind=c_size_t) :: sum_mem,free_mem,fourdp_smem
85 : : ! *************************************************************************
86 : :
87 : 0 : ider=0
88 [ # # ]: 0 : if (usepaw==1) then
89 : : nprojs=0
90 [ # # ]: 0 : do i = 1,natom
91 : 0 : nprojs = nprojs + cprj(i, 1)%nlmn
92 : : end do
93 : : end if
94 : :
95 : 0 : t_fft(1) = ngfft(3);
96 : 0 : t_fft(2) = ngfft(2);
97 : 0 : t_fft(3) = ngfft(1);
98 : :
99 : : #ifdef HAVE_GPU
100 : : call gpu_get_max_mem(free_mem)
101 : : free_mem = 0.85 * free_mem ! Cutting 15% out to be safe
102 : : #endif
103 : :
104 [ # # ]: 0 : do i=1,nband_k
105 [ # # ]: 0 : if(modulo(nband_k,i)/=0) cycle
106 : 0 : ndat_occ=nband_k/i
107 : 0 : sum_mem = 0
108 : :
109 : : ! cwavef
110 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*n4*n5*n6*ndat
111 [ # # ]: 0 : if(need_ghc) then
112 : : ! ghc
113 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*npw*ndat
114 : : ! ghc1
115 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*npw*ndat
116 : : ! ghc2
117 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*npw*ndat
118 : : end if
119 : :
120 : : ! cwavef_r
121 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*n4*n5*n6*ndat
122 : : ! vlocpsi_r
123 : 0 : sum_mem = sum_mem + INT(cplex_fock,c_size_t)*nfftf*ndat
124 : : ! work (ompgpu_fourwf internal array)
125 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*n4*n5*n6*ndat*ndat_occ
126 : : #ifdef HAVE_GPU
127 : : call gpu_fft_get_estimate_work_size(3, c_loc(t_fft), FFT_Z2Z, ndat*ndat_occ, fourdp_smem);
128 : : #endif
129 : 0 : sum_mem = sum_mem + fourdp_smem
130 : :
131 : : ! rhor_munu
132 : 0 : sum_mem = sum_mem + INT(cplex_fock,c_size_t)*nfftf*ndat_occ*ndat
133 : : ! rhog_munu
134 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nfftf*ndat_occ*ndat
135 : : ! vfock
136 : 0 : sum_mem = sum_mem + INT(cplex_fock,c_size_t)*nfftf*ndat_occ*ndat
137 : : ! occ
138 : 0 : sum_mem = sum_mem + INT(ndat_occ,c_size_t)
139 : :
140 : : !*Additional arrays in case of paw
141 [ # # ]: 0 : if (usepaw==1) then
142 [ # # # # ]: 0 : if ((optfor).and.(ieigen/=0)) then
143 : 0 : ider=3
144 : : ! forout
145 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*npw*nspinor*ndat_occ*ndat
146 : :
147 : : ! dprojs (gemm_nonlop_ompgpu internal work array)
148 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nprojs*npw*3
149 : : end if
150 : :
151 [ # # # # ]: 0 : if (optstr.and.(ieigen/=0)) then
152 : 0 : ider=3
153 : : ! strout
154 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*npw*nspinor*ndat_occ*ndat
155 : :
156 : : ! dprojs (gemm_nonlop_ompgpu internal work array)
157 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nprojs*npw*6
158 : : end if
159 : : ! grnhat_12
160 : 0 : ider=ider*2 ! Overestimate this buffer to ensure it fits as we don't manage OpenMP pool of GPU memory
161 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nfgd_max*nspinor**2*3*natom*(ider/3)*ndat_occ*ndat
162 : : ! gvnlxc
163 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*npw*nspinor*ndat_occ*ndat
164 : : ! rho12 (=nhat12 output of pawmknhat_psipsi_ndat)
165 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nfgd_max*nspinor**2*ndat_occ*ndat*natom
166 : : ! cprj1 (paw_psipsi internal work array)
167 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nprojs*nspinor*ndat
168 : : ! cprj2 (paw_psipsi internal work array)
169 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nprojs*nspinor*ndat
170 : : ! cpf (paw_psipsi internal work array)
171 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*lmn2_size*nspinor*ndat*ndat_occ*natom
172 : :
173 : : ! dijhat
174 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*lmn2_size*nspinor*ndat*ndat_occ*natom
175 : :
176 : : ! projs (gemm_nonlop_ompgpu internal work array)
177 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*nprojs*npw
178 : : end if
179 : :
180 : : ! cwaveocc_r
181 : 0 : sum_mem = sum_mem + INT(2,c_size_t)*n4*n5*n6*ndat_occ
182 : :
183 : 0 : sum_mem = sum_mem*dp
184 : :
185 [ # # ]: 0 : if(sum_mem < free_mem) exit
186 : : end do
187 : :
188 [ # # ]: 0 : if(sum_mem > free_mem) then
189 : 0 : ABI_WARNING("Test case doesn't fit in GPU memory. Try to lower bandpp.")
190 : : end if
191 : : !#ifdef DEBUG_VERBOSE
192 : 0 : write(std_out,*) "-----------DEBUG fock_getghc%select_ndat_occ_for_gpu : "
193 : 0 : write(std_out,'(A,F10.3,1x,A)') "Considered free GPU memory : ", real(free_mem,dp)/(1024*1024), "MiB"
194 : 0 : write(std_out,'(A,I4)') "selected ndat_occ : ", ndat_occ
195 : 0 : write(std_out,'(A,F10.3,1x,A)') "Forecasted consumed GPU memory : ", real(sum_mem,dp)/(1024*1024), "MiB"
196 : 0 : write(std_out,*) "-----------END DEBUG fock_getghc%select_ndat_occ_for_gpu : "
197 : : !#endif
198 : :
199 : 0 : end subroutine select_ndat_occ_for_gpu
200 : :
201 : : !!****f* ABINIT/fock_getghc
202 : : !! NAME
203 : : !! fock_getghc
204 : : !!
205 : : !! FUNCTION
206 : : !! Compute the matrix elements <G|Vx|psi> of the Fock operator.
207 : : !!
208 : : !! INPUTS
209 : : !! cwavef(2,npw*nspinor*ndat)= planewave coefficients of wavefunctions on which Fock operator is applied.
210 : : !! cwaveprj <type(pawcprj_type> = <cwavevf|proj>
211 : : !! gs_ham <type(gs_hamiltonian_type)>=all data for the Hamiltonian to be applied
212 : : !! mpi_enreg= information about MPI parallelization
213 : : !!
214 : : !! SIDE EFFECTS
215 : : !! ghc(2,npw*ndat)= matrix elements <G|H|C> or <G|H-lambda.S|C> (if sij_opt>=0 or =-1 in getghc)
216 : : !! contains the fock exchange term for cwavef at the end.
217 : : !!
218 : : !! NOTES
219 : : !! The current version assumes that:
220 : : !! * nspinor = 1
221 : : !! * no "my_nspinor"
222 : : !! * no restriction to the value of istwfk_bz (but must be tested in all case)
223 : : !! * all the data for the occupied states (cgocc_bz) are the same as those for the current states (cg)
224 : : !!
225 : : !! SOURCE
226 : :
227 [ - + - + : 29057 : subroutine fock_getghc(cwavef,cwaveprj,ghc,gs_ham,mpi_enreg,ndat)
- + ]
228 : :
229 : : !Arguments ------------------------------------
230 : : ! Scalars
231 : : integer,intent(in) :: ndat
232 : : type(MPI_type),intent(in) :: mpi_enreg
233 : : type(gs_hamiltonian_type),target,intent(inout) :: gs_ham
234 : : ! Arrays
235 : : type(pawcprj_type),intent(inout) :: cwaveprj(:,:)
236 : : real(dp),intent(inout) :: cwavef(:,:)!,ghc(2,gs_ham%npw_k)
237 : : real(dp),intent(inout) :: ghc(:,:)
238 : :
239 : : !Local variables-------------------------------
240 : : ! Scalars
241 : : integer,parameter :: tim_fourwf_fock_getghc=10,tim_fourdp_fock_getghc=10
242 : : integer :: bdtot_jindex,choice,cplex_fock,cplex_dij,cpopt,i1,i2,i3,ia,iatom,idat,idat_occ,idat_tot,iatm
243 : : integer :: iband_cprj,ider,idir,idir1,ier,ii,ind,ipw,ieigen,ifft,itypat,izero,jband,jbg,jcg,jkg
244 : : integer :: jkpt,my_jsppol,jstwfk,lmn2_size,mgfftf,mpw,n1,n2,n3,n4,n5,n6,ndat_occ,ndat_tot
245 : : integer :: n1f,n2f,n3f,n4f,n5f,n6f,natom,nband_k,ndij,nfft,nfftf,nfftotf,nhat12_grdim,nnlout
246 : : integer :: nfgd_max
247 : : integer :: npw,npwj,nspden_fock,nspinor,nkpg,paw_opt,signs,tim_nonlop,gpu_option
248 : : integer, save :: ncount=0
249 : : logical :: need_ghc,qeq0
250 : : real(dp),parameter :: weight1=one
251 : : real(dp) :: doti,eigen,imcwf,imcwocc,imvloc,invucvol,recwf,recwocc,revloc,wtk,esum,esumi
252 : : complex(dp) :: cinvucvol,cucvol
253 : : type(fock_common_type),pointer :: fockcommon
254 : : type(fock_BZ_type),pointer :: fockbz
255 : : ! Arrays
256 : : integer :: ngfft(18),ngfftf(18)
257 : 29057 : integer,pointer :: gboundf(:,:),kg_occ(:,:),gbound_kp(:,:)
258 : : real(dp) :: fockstr(6),qphon(3),qvec_j(3),tsec(2),gsc_dum(2,0),rhodum(2,1)
259 : : real(dp) :: rhodum0(0,1,1)
260 : 29057 : real(dp), allocatable :: dummytab(:,:),dijhat(:,:,:,:,:),dijhat_tmp(:,:,:),ffnl_kp_dum(:,:,:,:),kpg_kp(:,:),occ(:)
261 : 29057 : real(dp), allocatable, target :: gvnlxc(:,:),ghc1(:,:),ghc2(:,:),grnhat12(:,:,:,:,:,:),grnhat_12(:,:,:,:,:,:,:),forikpt(:,:,:)
262 : 29057 : real(dp), allocatable :: rho12(:,:,:,:,:,:),rhog_munu(:,:,:,:),rhor_munu(:,:,:,:),vlocpsi_r(:,:),strdat(:,:,:,:)
263 : 29057 : real(dp), allocatable :: vfock(:,:,:),psilocal(:,:,:),enlout_dum(:),vectin_dum(:,:),vqg(:),forout(:,:),strout(:,:),for1(:,:,:,:)
264 : 29057 : real(dp), allocatable,target ::cwavef_r(:,:,:,:),cwavef_rep(:,:),vdotr(:,:,:,:),vdoti(:),vfockstr(:,:,:)
265 : 29057 : real(dp), contiguous, pointer :: cwaveocc_r(:,:,:,:,:)
266 : 29057 : type(pawcprj_type),pointer :: cwaveocc_prj(:,:)
267 : 29057 : type(pawcprj_type),pointer :: cwaveocc_prj_rep(:,:)
268 : : real(dp) :: rprimd(3,3),for12(3)
269 : 29057 : integer, contiguous, pointer :: atom_ifftsph(:,:),atom_nfgd(:)
270 : 29057 : real(dp), contiguous, pointer :: stress_ikpt(:,:),atom_rfgd(:,:,:)
271 : : ! *************************************************************************
272 : :
273 : 29057 : ncount=ncount+1
274 : :
275 : 29057 : call timab(1504,1,tsec) ; call timab(1505,-1,tsec) ; call timab(1515,-1,tsec) ; call timab(1541,-1,tsec)
276 : :
277 [ - + ]: 29057 : ABI_CHECK(associated(gs_ham%fockcommon),"fock_common must be associated!")
278 : 29057 : fockcommon => gs_ham%fockcommon
279 [ - + ]: 29057 : ABI_CHECK(associated(gs_ham%fockbz),"fock_bz must be associated!")
280 : 29057 : fockbz => gs_ham%fockbz
281 : :
282 [ - + ]: 29057 : ABI_CHECK(gs_ham%nspinor==1,"only allowed for nspinor=1!")
283 [ - + ]: 29057 : ABI_CHECK(gs_ham%npw_k==gs_ham%npw_kp,"only allowed for npw_k=npw_kp (ground state)!")
284 [ + + ]: 29057 : if (fockcommon%usepaw==1) then
285 [ + - - + ]: 2115 : ABI_CHECK((size(cwaveprj,1)==gs_ham%natom.and.size(cwaveprj,2)==gs_ham%nspinor*ndat),"error on cwaveprj dims")
286 : : end if
287 : 29057 : need_ghc=(size(ghc,2)>0)
288 : :
289 : : !Some constants
290 : 29057 : invucvol=1.d0/sqrt(gs_ham%ucvol)
291 : 29057 : cinvucvol=dcmplx(invucvol,0.0_dp)
292 : 29057 : cucvol=dcmplx(sqrt(gs_ham%ucvol),0.0_dp)
293 : 29057 : call matr3inv(gs_ham%gprimd,rprimd)
294 : 29057 : cplex_fock=2;nspden_fock=1
295 : 29057 : natom=fockcommon%natom
296 : 29057 : nspinor=gs_ham%nspinor
297 : : mpw=maxval(fockbz%npwarr)
298 : 29057 : npw=gs_ham%npw_k
299 : 29057 : gpu_option=gs_ham%gpu_option
300 : 29057 : ider=0;izero=0
301 [ + + ]: 29057 : if (fockcommon%usepaw==1) then
302 [ + + ]: 40185 : nfft =fockcommon%pawfgr%nfftc ; ngfft =fockcommon%pawfgr%ngfftc
303 [ + + ]: 40185 : nfftf=fockcommon%pawfgr%nfft ; ngfftf=fockcommon%pawfgr%ngfft
304 : 2115 : mgfftf=fockcommon%pawfgr%mgfft
305 : : else
306 : 26942 : nfft =gs_ham%nfft ; nfftf =nfft
307 [ + + ]: 511898 : ngfft=gs_ham%ngfft ; ngfftf=ngfft
308 : 26942 : mgfftf=gs_ham%mgfft
309 : : end if
310 : 29057 : n1=ngfft(1);n2=ngfft(2);n3=ngfft(3)
311 : 29057 : n4=ngfft(4);n5=ngfft(5);n6=ngfft(6)
312 : 29057 : n1f=ngfftf(1);n2f=ngfftf(2);n3f=ngfftf(3)
313 : 29057 : n4f=ngfftf(4);n5f=ngfftf(5);n6f=ngfftf(6)
314 : :
315 : : !*Max number of fine-grid points in a PAW augmentation sphere (over all atoms).
316 : : !*Used to size the small per-atom rho12/grnhat_12 sphere buffers instead of the full FFT grid.
317 : 29057 : nfgd_max=1
318 [ + + ]: 29057 : if (fockcommon%usepaw==1) then
319 [ + + ]: 6345 : nfgd_max=max(1,maxval(fockcommon%pawfgrtab(1:natom)%nfgd))
320 : : end if
321 : :
322 : : ! ===========================
323 : : ! === Initialize arrays ===
324 : : ! ===========================
325 : : ! transient optfor and optstress
326 : : ! fockcommon%optfor=.false.
327 : : ! fockcommon%optstr=.false.
328 : : !*Initialization of local pointers
329 : : !*Initialization of the array cwavef_r
330 : : !*cwavef_r = current wavefunction in r-space
331 [ - + + - : 145285 : ABI_MALLOC(cwavef_r,(2,n4f,n5f,n6f*ndat))
- + + - -
+ + - + -
+ - - + -
+ - + ]
332 : : #ifdef HAVE_OPENMP_OFFLOAD
333 : : !$OMP TARGET ENTER DATA MAP(alloc:cwavef_r) IF(gpu_option==ABI_GPU_OPENMP)
334 : : #endif
335 : : !*dummytab = variables for fourwf
336 [ - + + - : 87171 : ABI_MALLOC(dummytab,(2,nfft*ndat))
+ - + - -
+ - + ]
337 : : !*vqg = 4pi/(G+q)**2
338 [ - + + - : 87171 : ABI_MALLOC(vqg,(nfftf))
+ - + - -
+ - + ]
339 : :
340 [ + + ]: 29057 : if(need_ghc) then
341 : : !*Initialization of the array ghc1
342 : : !*ghc1 will contain the exact exchange contribution to the Hamiltonian
343 [ - + + - : 80943 : ABI_MALLOC(ghc1,(2,npw*ndat))
+ - + - -
+ - + -
+ ]
344 [ - + + - : 53962 : ABI_MALLOC(ghc2,(2,npw*ndat))
- + - + -
+ ]
345 : : #ifdef HAVE_OPENMP_OFFLOAD
346 : : !$OMP TARGET ENTER DATA MAP(alloc:ghc1,ghc2,vqg) IF(gpu_option==ABI_GPU_OPENMP)
347 : : #endif
348 [ + - ]: 26981 : if(gpu_option==ABI_GPU_DISABLED) then
349 [ + + + + ]: 20059469 : ghc1=zero
350 [ + + + + ]: 20059469 : ghc2=zero
351 [ # # ]: 0 : else if(gpu_option==ABI_GPU_OPENMP) then
352 : 0 : call gpu_set_to_zero(ghc1, int(2,c_size_t)*npw*ndat)
353 : 0 : call gpu_set_to_zero(ghc2, int(2,c_size_t)*npw*ndat)
354 : : end if
355 : : end if
356 : : !*Initialization of the array vlocpsi_r
357 : : !*vlocpsi_r = partial local Fock operator applied to cwavef in r-space and summed over all occupied (jkpt,mu)
358 [ - + + - : 116228 : ABI_MALLOC(vlocpsi_r,(cplex_fock*nfftf,ndat))
- + + - +
- + - - +
- + ]
359 : : #ifdef HAVE_OPENMP_OFFLOAD
360 : : !$OMP TARGET ENTER DATA MAP(alloc:vlocpsi_r) IF(gpu_option==ABI_GPU_OPENMP)
361 : : #endif
362 [ + - ]: 29057 : if(gpu_option==ABI_GPU_DISABLED) then
363 [ + + + + ]: 457400022 : vlocpsi_r=zero
364 [ # # ]: 0 : else if(gpu_option==ABI_GPU_OPENMP) then
365 : 0 : call gpu_set_to_zero(vlocpsi_r, int(cplex_fock,c_size_t)*nfftf*ndat)
366 : : end if
367 : :
368 : : !*Additional arrays in case of paw
369 [ + + ]: 29057 : if (fockcommon%usepaw==1) then
370 : 2115 : nhat12_grdim=0
371 : : end if
372 : :
373 [ + + ]: 29057 : if (fockcommon%usepaw==1) then
374 [ + + + - ]: 2115 : if ((fockcommon%optfor).and.(fockcommon%ieigen/=0)) then
375 [ - + + - : 4304 : ABI_MALLOC(forikpt,(3,natom,ndat))
- + + - +
- + - - +
- + ]
376 [ + + + + : 10760 : forikpt=zero
+ + ]
377 : : end if
378 : : end if
379 [ + + + + ]: 29057 : if (fockcommon%usepaw==1.or.fockcommon%optstr) then
380 [ + - - + : 15812 : ABI_MALLOC(gboundf,(2*mgfftf+8,2))
+ - + - +
- + - - +
- + ]
381 [ + - - - : 3953 : call sphereboundary(gboundf,gs_ham%istwf_k,gs_ham%kg_k,mgfftf,npw)
- - - - -
+ - + ]
382 : : else
383 : 25104 : gboundf=>gs_ham%gbound_k
384 : : end if
385 : : ! ==========================================
386 : : ! === Get cwavef in real space using FFT ===
387 : : ! ==========================================
388 [ + - ]: 29057 : if(gpu_option==ABI_GPU_DISABLED) then
389 [ + + + + : 741148478 : cwavef_r=zero
+ + + + ]
390 [ # # ]: 0 : else if(gpu_option==ABI_GPU_OPENMP) then
391 : 0 : call gpu_set_to_zero(cwavef_r, int(2,c_size_t)*n4f*n5f*n6f*ndat)
392 : : end if
393 : 29057 : call timab(1515,2,tsec) ; call timab(1541,-2,tsec) ; call timab(1512,-1,tsec)
394 : : call fourwf(0,rhodum0,cwavef,rhodum,cwavef_r,gboundf,gboundf,gs_ham%istwf_k,gs_ham%kg_k,gs_ham%kg_k,&
395 : : & mgfftf,mpi_enreg,ndat,ngfftf,npw,1,n4f,n5f,n6f,0,tim_fourwf_fock_getghc,weight1,weight1,&
396 [ - + - - : 29057 : & gpu_option=gs_ham%gpu_option)
- - - - +
- - - - -
- - + - -
- - - - -
- + - - -
- - + - +
- + - + ]
397 : 29057 : call timab(1512,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1541,-1,tsec)
398 [ + - ]: 29057 : if(gpu_option==ABI_GPU_DISABLED) then
399 [ + + + + : 741148478 : cwavef_r=cwavef_r*invucvol
+ + + + ]
400 [ # # ]: 0 : else if(gpu_option==ABI_GPU_OPENMP) then
401 : 0 : call abi_xscal(n4f*n5f*n6f*ndat,cinvucvol,cwavef_r,1,x_cplx=2,gpu_option=gpu_option)
402 : : end if
403 : :
404 : : ! =====================================================
405 : : ! === Select the states in cgocc_bz with the same spin ===
406 : : ! =====================================================
407 : : !* Initialization of the indices/shifts, according to the value of isppol
408 : : !* bdtot_jindex = shift to be applied on the location of data in the array occ_bz ?
409 : 29057 : bdtot_jindex=0
410 : : !* jbg = shift to be applied on the location of data in the array cprj/occ
411 : 29057 : jbg=0;jcg=0
412 : 29057 : my_jsppol=fockcommon%isppol
413 [ + + - + ]: 29057 : if((fockcommon%isppol==2).and.(mpi_enreg%nproc_spkpt/=1)) my_jsppol=1
414 : :
415 : : !===================================
416 : : !=== Loop on the k-points in IBZ ===
417 : : !===================================
418 : 29057 : jkg=0
419 : :
420 [ + - ]: 29057 : if (associated(gs_ham%ph3d_kp)) then
421 : 29057 : nullify (gs_ham%ph3d_kp)
422 : : end if
423 : :
424 : 29057 : call timab(1505,2,tsec) ; call timab(1506,-1,tsec) ; call timab(1541,-2,tsec)
425 : :
426 [ + + ]: 232015 : do jkpt=1,fockbz%mkpt
427 : :
428 [ - + ]: 202958 : if(fockbz%nbandocc_bz(jkpt,my_jsppol)==0) cycle
429 : :
430 : 202958 : call timab(1521,1,tsec)
431 : :
432 : : !* nband_k = number of bands at point k_j
433 : 202958 : nband_k=fockbz%nbandocc_bz(jkpt,my_jsppol)
434 : :
435 : : ! Select ndat_occ value :
436 : : ! GPU : check GPU memory available and compute maximum value for ndat_occ
437 : : ! CPU : 4, because it seems to be optimal from my observations
438 [ - + ]: 202958 : if(gpu_option/=ABI_GPU_DISABLED) then
439 : 0 : lmn2_size=0
440 [ # # ]: 0 : if(fockcommon%usepaw==1) lmn2_size=fockcommon%pawtab(1)%lmn2_size
441 : : call select_ndat_occ_for_gpu(ndat_occ,nband_k,ndat,npw,cplex_fock,&
442 : : & nfftf,ngfftf,n4f,n5f,n6f,natom,nspinor,lmn2_size,&
443 [ # # # # : 0 : & fockcommon%usepaw,cwaveprj,fockcommon%ieigen,need_ghc,fockcommon%optfor,fockcommon%optstr,nfgd_max)
# # # # #
# ]
444 : : else
445 : 202958 : ndat_occ=min(nband_k,16)
446 [ + - ]: 202958 : do ii=1,nband_k
447 [ + - ]: 202958 : if(modulo(nband_k,ndat_occ)==0) exit
448 : 202958 : ndat_occ=ndat_occ-1
449 : : end do
450 : : end if
451 : : !* ndat_tot = combined batch size (ndat*ndat_occ), used to batch nonlop calls
452 : 202958 : ndat_tot=ndat*ndat_occ
453 : :
454 : : !* wtk = weight in BZ of this k point
455 : 202958 : wtk=fockbz%wtk_bz(jkpt) !*sqrt(gs_ham%ucvol)
456 : : !* jstwfk= how is stored the wavefunction
457 : 202958 : jstwfk=fockbz%istwfk_bz(jkpt)
458 : : !* npwj= number of plane wave in basis for the wavefunction
459 : 202958 : npwj=fockbz%npwarr(jkpt)
460 : : !* Basis sphere of G vectors
461 [ + - ]: 202958 : if (allocated(fockbz%cgocc)) then
462 : 202958 : gbound_kp => fockbz%gbound_bz(:,:,jkpt)
463 : 202958 : kg_occ => fockbz%kg_bz(:,1+jkg:npwj+jkg)
464 : : end if
465 : :
466 : : !*rhormunu = overlap matrix between cwavef and (jkpt,mu) in R-space
467 [ - + + - : 1217748 : ABI_MALLOC(rhor_munu,(cplex_fock,nfftf,ndat_occ,ndat))
- + + - -
+ + - - +
+ - + - +
- - + -
+ ]
468 : : #ifdef HAVE_OPENMP_OFFLOAD
469 : : !$OMP TARGET ENTER DATA MAP(alloc:rhor_munu) IF(gpu_option==ABI_GPU_OPENMP)
470 : : #endif
471 : : !*rhogmunu = overlap matrix between cwavef and (jkpt,mu) in G-space
472 [ - + + - : 1014790 : ABI_MALLOC(rhog_munu,(2,nfftf,ndat_occ,ndat))
- + + - -
+ + - + -
+ - - + -
+ ]
473 : : #ifdef HAVE_OPENMP_OFFLOAD
474 : : !$OMP TARGET ENTER DATA MAP(alloc:rhog_munu) IF(gpu_option==ABI_GPU_OPENMP)
475 : : #endif
476 : : !*vfock = Fock potential
477 [ - + + - : 1014790 : ABI_MALLOC(vfock,(cplex_fock*nfftf,ndat_occ,ndat))
- + + - -
+ + - + -
+ - - + -
+ ]
478 : : #ifdef HAVE_OPENMP_OFFLOAD
479 : : !$OMP TARGET ENTER DATA MAP(alloc:vfock) IF(gpu_option==ABI_GPU_OPENMP)
480 : : #endif
481 [ - + + - : 608874 : ABI_MALLOC(occ,(ndat_occ))
+ - + - -
+ - + ]
482 : : #ifdef HAVE_OPENMP_OFFLOAD
483 : : !$OMP TARGET ENTER DATA MAP(alloc:occ) IF(gpu_option==ABI_GPU_OPENMP)
484 : : #endif
485 : :
486 : : !*Additional arrays in case of paw
487 [ + + ]: 202958 : if (fockcommon%usepaw==1) then
488 [ + + + - ]: 10110 : if ((fockcommon%optfor).and.(fockcommon%ieigen/=0)) then
489 : 2996 : ider=3
490 [ - + + - : 8988 : ABI_MALLOC(forout,(2,npw*nspinor*ndat_tot))
+ - + - -
+ - + -
+ ]
491 : : #ifdef HAVE_OPENMP_OFFLOAD
492 : : !$OMP TARGET ENTER DATA MAP(alloc:forout) IF(gpu_option==ABI_GPU_OPENMP)
493 : : #endif
494 : : end if
495 : :
496 [ + + + - ]: 10110 : if (fockcommon%optstr.and.(fockcommon%ieigen/=0)) then
497 : 1870 : ider=3
498 [ - + + - : 5610 : ABI_MALLOC(strout,(2,npw*nspinor*ndat_tot))
+ - + - -
+ - + -
+ ]
499 : : #ifdef HAVE_OPENMP_OFFLOAD
500 : : !$OMP TARGET ENTER DATA MAP(alloc:strout) IF(gpu_option==ABI_GPU_OPENMP)
501 : : #endif
502 : : end if
503 [ - + + - : 80880 : ABI_MALLOC(grnhat_12,(2,nfgd_max,nspinor**2,3,natom*(ider/3),ndat_occ,ndat))
- + + - +
- + + + -
- + + - -
+ + - + -
+ + - + -
+ - + ]
504 [ - + + - : 30330 : ABI_MALLOC(gvnlxc,(2,npw*nspinor*ndat_tot))
+ - + - -
+ - + -
+ ]
505 [ - + - + ]: 10110 : ABI_MALLOC(grnhat12,(1,1,1,1,1,1))
506 [ - + + - : 70770 : ABI_MALLOC(rho12,(2,nfgd_max,nspinor**2,ndat_occ,ndat,natom))
- + + - -
+ + - - +
+ - - + +
- + - + -
- + - + -
+ ]
507 : :
508 [ - + + - : 30330 : ABI_MALLOC(atom_nfgd, (natom))
+ - + - -
+ - + ]
509 [ + + ]: 30330 : do iatom=1,natom
510 : 30330 : atom_nfgd(iatom) = fockcommon%pawfgrtab(iatom)%nfgd
511 : : end do
512 [ - + + - : 40440 : ABI_MALLOC(atom_ifftsph, (nfgd_max, natom))
- + + - +
- + - - +
- + ]
513 [ - + + - : 40440 : ABI_MALLOC(atom_rfgd, (3, nfgd_max, natom))
- + + - +
- + - - +
- + ]
514 [ + + ]: 30330 : do iatom=1,natom
515 [ + + ]: 16123784 : atom_ifftsph(1:atom_nfgd(iatom),iatom) = fockcommon%pawfgrtab(iatom)%ifftsph(1:atom_nfgd(iatom))
516 [ + + + + ]: 64444586 : atom_rfgd(:,1:atom_nfgd(iatom),iatom) = fockcommon%pawfgrtab(iatom)%rfgd(:,1:atom_nfgd(iatom))
517 : : end do
518 : :
519 : : #ifdef HAVE_OPENMP_OFFLOAD
520 : : !$OMP TARGET ENTER DATA MAP(alloc:grnhat_12) IF(gpu_option==ABI_GPU_OPENMP .and. ider==3)
521 : : !$OMP TARGET ENTER DATA MAP(alloc:gvnlxc) IF(gpu_option==ABI_GPU_OPENMP)
522 : : !$OMP TARGET ENTER DATA MAP(alloc:rho12) IF(gpu_option==ABI_GPU_OPENMP)
523 : : !$OMP TARGET ENTER DATA MAP(alloc:atom_ifftsph,atom_nfgd,atom_rfgd) IF(gpu_option==ABI_GPU_OPENMP)
524 : : !$OMP TARGET UPDATE TO(atom_ifftsph,atom_nfgd,atom_rfgd) IF(gpu_option==ABI_GPU_OPENMP)
525 : : #endif
526 : : end if
527 : :
528 : : !* Load k^prime hamiltonian in the gs_ham datastructure
529 : : ! Note: ffnl_kp / ph3d_kp / gbound_kp are not used
530 : :
531 [ + - ]: 202958 : if (.not. associated(gs_ham%ph3d_kp)) then
532 [ - + + - : 811832 : ABI_MALLOC(gs_ham%ph3d_kp,(2,npwj,gs_ham%matblk))
- + + - +
- + - - +
- + ]
533 : : end if
534 : :
535 : : call gs_ham%load_kprime(kpt_kp=fockbz%kptns_bz(:,jkpt),&
536 : 202958 : & istwf_kp=jstwfk,npw_kp=npwj,kg_kp=fockbz%kg_bz(:,1+jkg:npwj+jkg))
537 : : !* Some temporary allocations needed for PAW
538 [ + + ]: 202958 : if (fockcommon%usepaw==1) then
539 [ - + + - : 30330 : ABI_MALLOC(enlout_dum,(ndat_tot))
+ - + - -
+ - + -
+ ]
540 [ - + + - : 30330 : ABI_MALLOC(vectin_dum,(2,npwj*nspinor*ndat_tot))
+ - + - -
+ - + -
+ ]
541 [ + + + + ]: 29090475 : vectin_dum=zero
542 [ - + + - : 50550 : ABI_MALLOC(ffnl_kp_dum,(npwj,1,gs_ham%lmnmax,gs_ham%ntypat))
- + + - -
+ + - + -
+ - - + -
+ - + ]
543 : 10110 : nkpg=size(gs_ham%kpg_k,2)
544 [ - + + - : 40440 : ABI_MALLOC(kpg_kp,(npwj,nkpg))
- + + - +
- + - - +
- + - + ]
545 [ + - ]: 10110 : if (nkpg>0) then
546 [ + - - - : 10110 : call mkkpg(gs_ham%kg_kp,kpg_kp,gs_ham%kpt_kp,nkpg,npwj)
- - - - -
+ ]
547 : : end if
548 : 10110 : call gs_ham%load_kprime(ffnl_kp=ffnl_kp_dum,kpg_kp=kpg_kp)
549 : : end if
550 : :
551 : : ! ======================================
552 : : ! === Calculate the vector q=k_i-k_j ===
553 : : ! ======================================
554 : : !* Evaluation of kpoint_j, the considered k-point in reduced coordinates
555 : : ! kpoint_j(:)=fockbz%kptns_bz(:,jkpt)
556 : : !* the vector qvec is expressed in reduced coordinates.
557 : : ! qvec(:)=kpoint_i(:)-kpoint_j(:)
558 [ + + ]: 811832 : qvec_j(:)=gs_ham%kpt_k(:)-fockbz%kptns_bz(:,jkpt)
559 : 202958 : qeq0=(qvec_j(1)**2+qvec_j(2)**2+qvec_j(3)**2<1.d-15)
560 : :
561 : : ! Get the Coulomb interaction in reciprocal space
562 : 202958 : call bare_vqg(qvec_j,fockcommon,gs_ham%gmet,nfftf,fockbz%nkpt_bz,ngfftf,gs_ham%ucvol,vqg)
563 : : #ifdef HAVE_OPENMP_OFFLOAD
564 : : !$OMP TARGET UPDATE TO(vqg) IF(gpu_option==ABI_GPU_OPENMP)
565 : : #endif
566 : :
567 : 202958 : call timab(1521,2,tsec)
568 : :
569 : : ! =================================================
570 : : ! === Loop on the band indices jband of cgocc_k ===
571 : : ! =================================================
572 [ - + - - : 405916 : do jband=1,nband_k,ndat_occ
- + - + ]
573 : :
574 : : !* occ = occupancy of jband at this k point
575 [ + + ]: 1117532 : occ(1:ndat_occ)=fockbz%occ_bz(jband+bdtot_jindex:jband+ndat_occ-1+bdtot_jindex,my_jsppol)
576 [ + - - + : 1320490 : if(maxval(occ(1:ndat_occ))<tol8) cycle
- - + + -
+ ]
577 : : #ifdef HAVE_OPENMP_OFFLOAD
578 : : !$OMP TARGET UPDATE TO(occ) IF(gpu_option==ABI_GPU_OPENMP)
579 : : #endif
580 : :
581 : : ! This timing is placed after the cycle ...
582 : 202958 : call timab(1522,1,tsec) ; call timab(1542,-1,tsec)
583 : :
584 : : ! ==============================================
585 : : ! === Get cwaveocc_r in real space using FFT ===
586 : : ! ==============================================
587 [ - + ]: 202958 : if (allocated(fockbz%cwaveocc_bz)) then
588 : 0 : cwaveocc_r => fockbz%cwaveocc_bz(:,:,:,:,jband+jbg:jband+jbg+ndat_occ-1,my_jsppol)
589 : : else
590 [ - + + - : 1217748 : ABI_MALLOC(cwaveocc_r,(2,n4f,n5f,n6f,ndat_occ))
- + + - -
+ + - - +
+ - + - +
- - + -
+ ]
591 : : #ifdef HAVE_OPENMP_OFFLOAD
592 : : !$OMP TARGET ENTER DATA MAP(alloc:cwaveocc_r) IF(gpu_option==ABI_GPU_OPENMP)
593 : : #endif
594 [ + - ]: 202958 : if(gpu_option==ABI_GPU_DISABLED) then
595 [ + + + + : 13580932253 : cwaveocc_r=zero
+ + + + +
+ ]
596 [ # # ]: 0 : else if(gpu_option==ABI_GPU_OPENMP) then
597 : 0 : call gpu_set_to_zero(cwaveocc_r, int(2,c_size_t)*n4f*n5f*n6f*ndat_occ)
598 : : end if
599 : 202958 : call timab(1515,2,tsec) ; call timab(1512,-1,tsec) ; call timab(1542,-2,tsec)
600 : : call fourwf(1,rhodum0,fockbz%cgocc(:,1+jcg+npwj*(jband-1):jcg+(jband+ndat_occ-1)*npwj,my_jsppol),rhodum,cwaveocc_r, &
601 : : & gbound_kp,gbound_kp,jstwfk,kg_occ,kg_occ,mgfftf,mpi_enreg,ndat_occ,ngfftf,&
602 [ - + - + : 202958 : & npwj,1,n4f,n5f,n6f,0,tim_fourwf_fock_getghc,weight1,weight1,gpu_option=gs_ham%gpu_option)
- + - + ]
603 : 202958 : call timab(1512,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1542,-1,tsec)
604 [ + - ]: 202958 : if(gpu_option==ABI_GPU_DISABLED) then
605 [ + + + + : 13580932253 : cwaveocc_r=cwaveocc_r*invucvol
+ + + + +
+ ]
606 [ # # ]: 0 : else if(gpu_option==ABI_GPU_OPENMP) then
607 : 0 : call abi_xscal(n4f*n5f*n6f*ndat_occ,cinvucvol,cwaveocc_r,1,x_cplx=2,gpu_option=gpu_option)
608 : : end if
609 : : end if
610 : :
611 : : ! ================================================
612 : : ! === Get the overlap density matrix rhor_munu ===
613 : : ! ================================================
614 : : !* Calculate the overlap density matrix in real space = conj(cwaveocc_r)*cwavef_r
615 : : !* rhor_munu will contain the overlap density matrix.
616 : : ! vfock=-int{conj(cwaveocc_r)*cwavef_r*dr'/|r-r'|}
617 : :
618 : 202958 : call timab(1522,2,tsec) ; call timab(1542,-2,tsec) ; call timab(1523,-1,tsec)
619 : :
620 [ + - ]: 202958 : if(gpu_option==ABI_GPU_DISABLED) then
621 : : !$OMP PARALLEL DO COLLAPSE(2) &
622 : : !$OMP& PRIVATE(ind,imcwf,recwf,recwocc,imcwocc)
623 [ + + ]: 405916 : do idat=1,ndat
624 [ + + ]: 1320490 : do idat_occ=1,ndat_occ
625 [ + + ]: 15610575 : do i3=1,n3f
626 [ + + ]: 251580668 : do i2=1,n2f
627 [ + + ]: 4356612612 : do i1=1,n1f
628 : 4105946518 : ind=i1+(i2-1)*n1f+(i3-1)*n2f*n1f
629 : 4105946518 : recwf =cwavef_r(1,i1,i2,(idat-1)*n3f+i3)
630 : 4105946518 : imcwf =cwavef_r(2,i1,i2,(idat-1)*n3f+i3)
631 : 4105946518 : recwocc=cwaveocc_r(1,i1,i2,i3,idat_occ)
632 : 4105946518 : imcwocc=cwaveocc_r(2,i1,i2,i3,idat_occ)
633 : 4105946518 : rhor_munu(1,ind,idat_occ,idat)= recwocc*recwf+imcwocc*imcwf
634 : 4342119569 : rhor_munu(2,ind,idat_occ,idat)= recwocc*imcwf-imcwocc*recwf
635 : : end do ! i1
636 : : end do ! i2
637 : : end do ! i3
638 : : end do ! idat_occ
639 : : end do ! idat
640 : : else if(gpu_option==ABI_GPU_OPENMP) then
641 : : #ifdef HAVE_OPENMP_OFFLOAD
642 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(2) &
643 : : !$OMP& MAP(to:rhor_munu,cwavef_r,cwaveocc_r) PRIVATE(idat,idat_occ)
644 : : do idat=1,ndat
645 : : do idat_occ=1,ndat_occ
646 : : !$OMP PARALLEL DO COLLAPSE(3) PRIVATE(ind,i3,i2,i1)
647 : : do i3=1,n3f
648 : : do i2=1,n2f
649 : : do i1=1,n1f
650 : : ind=i1+(i2-1)*n1f+(i3-1)*n2f*n1f
651 : : rhor_munu(1,ind,idat_occ,idat)= cwaveocc_r(1,i1,i2,i3,idat_occ)*cwavef_r(1,i1,i2,(idat-1)*n3f+i3)&
652 : : +cwaveocc_r(2,i1,i2,i3,idat_occ)*cwavef_r(2,i1,i2,(idat-1)*n3f+i3)
653 : : rhor_munu(2,ind,idat_occ,idat)= cwaveocc_r(1,i1,i2,i3,idat_occ)*cwavef_r(2,i1,i2,(idat-1)*n3f+i3)&
654 : : -cwaveocc_r(2,i1,i2,i3,idat_occ)*cwavef_r(1,i1,i2,(idat-1)*n3f+i3)
655 : : end do ! i1
656 : : end do ! i2
657 : : end do ! i3
658 : : end do ! idat_occ
659 : : end do ! idat
660 : : #endif
661 : : end if ! gpu_option
662 : :
663 : 202958 : call timab(1523,2,tsec)
664 : :
665 : : ! =======================================================
666 : : ! === Add compensation charge density in the PAW case ===
667 : : ! =======================================================
668 : :
669 : 202958 : call timab(1524,-1,tsec) ; call timab(1544,-1,tsec)
670 : :
671 [ + + ]: 202958 : if (fockcommon%usepaw==1) then
672 : :
673 : 10110 : iband_cprj=(my_jsppol-1)*fockbz%mkptband+jbg+jband
674 : 10110 : cwaveocc_prj=>fockbz%cwaveocc_prj(:,iband_cprj:iband_cprj+ndat_occ*nspinor-1)
675 : :
676 : : call pawmknhat_psipsi_ndat(cwaveprj(:,:),cwaveocc_prj(:,:),&
677 : : & ider,izero,natom,natom,nfftf,ngfftf,&
678 : : & nhat12_grdim,nspinor,fockcommon%ntypat,ndat,ndat_occ,fockbz%pawang,fockcommon%pawfgrtab,grnhat12,&
679 : : & rho12,gs_ham%nattyp,&
680 : : & fockcommon%pawtab,gprimd=gs_ham%gprimd,grnhat_12=grnhat_12,qphon=qvec_j,&
681 [ - + - - : 10110 : & xred=gs_ham%xred,atindx=gs_ham%atindx,gpu_option=gpu_option)
- - - - +
- - - - -
- - - + -
+ - + -
+ ]
682 : :
683 [ + - ]: 10110 : if(gpu_option==ABI_GPU_DISABLED) then
684 : : !$OMP PARALLEL DO COLLAPSE(2) &
685 : : !$OMP& PRIVATE(idat,idat_occ,iatom,ifft,ind)
686 [ + + ]: 20220 : do idat=1,ndat
687 [ + + ]: 60321 : do idat_occ=1,ndat_occ
688 [ + + ]: 130413 : do iatom=1,natom
689 [ + + ]: 63588659 : do ifft=1,atom_nfgd(iatom)
690 : 63468356 : ind = atom_ifftsph(ifft,iatom)
691 : 63468356 : rhor_munu(1,ind,idat_occ,idat) = rhor_munu(1,ind,idat_occ,idat) + rho12(1,ifft,nspinor,idat_occ,idat,iatom)
692 : 63548558 : rhor_munu(2,ind,idat_occ,idat) = rhor_munu(2,ind,idat_occ,idat) - rho12(2,ifft,nspinor,idat_occ,idat,iatom)
693 : : end do ! ifft
694 : : end do ! iatom
695 : : end do ! idat_occ
696 : : end do ! idat
697 : : else if(gpu_option==ABI_GPU_OPENMP) then
698 : : #ifdef HAVE_OPENMP_OFFLOAD
699 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(2) &
700 : : !$OMP& MAP(to:rhor_munu,rho12,atom_nfgd,atom_ifftsph) PRIVATE(idat,idat_occ)
701 : : do idat=1,ndat
702 : : do idat_occ=1,ndat_occ
703 : : !$OMP PARALLEL DO PRIVATE(iatom,ifft,ind)
704 : : do iatom=1,natom
705 : : do ifft=1,atom_nfgd(iatom)
706 : : ind = atom_ifftsph(ifft,iatom)
707 : : rhor_munu(1,ind,idat_occ,idat) = rhor_munu(1,ind,idat_occ,idat) + rho12(1,ifft,nspinor,idat_occ,idat,iatom)
708 : : rhor_munu(2,ind,idat_occ,idat) = rhor_munu(2,ind,idat_occ,idat) - rho12(2,ifft,nspinor,idat_occ,idat,iatom)
709 : : end do ! ifft
710 : : end do ! iatom
711 : : end do ! idat_occ
712 : : end do ! idat
713 : : #endif
714 : : end if
715 : : end if
716 : :
717 : :
718 : 202958 : call timab(1515,2,tsec) ; call timab(1513,-1,tsec) ; call timab(1544,-2,tsec)
719 : : ! Perform an FFT using fourwf to get rhog_munu = FFT^-1(rhor_munu)
720 : : call fourdp(cplex_fock,rhog_munu,rhor_munu,-1,mpi_enreg,nfftf,ndat*ndat_occ,&
721 : 202958 : & ngfftf,tim_fourdp_fock_getghc,gpu_option=gpu_option)
722 : 202958 : call timab(1513,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1544,-1,tsec)
723 : :
724 [ + + + - ]: 202958 : if(fockcommon%optstr.and.(fockcommon%ieigen/=0)) then
725 [ - + + - : 51840 : ABI_MALLOC(vfockstr, (6,ndat_occ,ndat))
- + + - +
- + - - +
- + ]
726 : : call strfock(fockcommon,gs_ham%gprimd,vfockstr,&
727 : : & mpi_enreg,nfftf,ngfftf,fockbz%nkpt_bz,ndat*ndat_occ,rhog_munu,gs_ham%ucvol,&
728 : 12960 : & qvec_j,gpu_option=gpu_option)
729 [ + + ]: 25920 : do idat=1,ndat
730 [ + + ]: 87080 : do idat_occ=1,ndat_occ
731 [ + - + + : 808040 : fockcommon%stress_ikpt(:,fockcommon%ieigen+idat-1)=fockcommon%stress_ikpt(:,fockcommon%ieigen+idat-1)+vfockstr(:,idat_occ,idat)*occ(idat_occ)*wtk
+ + ]
732 : : end do ! idat_occ
733 : : end do ! idat
734 [ - + ]: 12960 : ABI_FREE(vfockstr)
735 [ + + + - ]: 12960 : if (fockcommon%usepaw==0.and.(.not.need_ghc)) then
736 [ + - ]: 11090 : if (allocated(fockbz%cgocc)) then
737 : : #ifdef HAVE_OPENMP_OFFLOAD
738 : : !$OMP TARGET EXIT DATA MAP(delete:cwaveocc_r) IF(gpu_option==ABI_GPU_OPENMP)
739 : : #endif
740 [ - + ]: 11090 : ABI_FREE(cwaveocc_r)
741 : : end if
742 : 11090 : call timab(1524,2,tsec) ; call timab(1544,-2,tsec)
743 : 11090 : cycle
744 : : end if
745 : : end if
746 : 191868 : call timab(1524,2,tsec) ; call timab(1544,-2,tsec)
747 : :
748 : : ! ===================================================
749 : : ! === Calculate the local potential vfockloc_munu ===
750 : : ! ===================================================
751 : : !* Apply the Poisson solver to "rhog_munu" while taking into account the effect of the vector "qvec"
752 : : !* This is precisely what is done in the subroutine hartre, with option cplex=2.
753 : : !* vfock will contain the local Fock potential, the result of hartre routine.
754 : : !* vfock = FFT( rhog_munu/|g+qvec|^2 )
755 : 191868 : call timab(1525,-1,tsec) ; call timab(1545,-1,tsec)
756 : : #if 0
757 : :
758 : : do idat=1,ndat
759 : : do idat_occ=1,ndat_occ
760 : : call timab(1515,-2,tsec) ; call timab(1513,-1,tsec)
761 : : call hartre(cplex_fock,fockcommon%gsqcut,fockcommon%usepaw,mpi_enreg,nfftf,ngfftf,&
762 : : & mpi_enreg%paral_kgb,rhog_munu(:,:,idat_occ,idat),rprimd,vfock(:,idat_occ,idat),divgq0=fock%divgq0,qpt=qvec_j)
763 : : call timab(1513,2,tsec) ; call timab(1515,-1,tsec)
764 : : end do ! idat_occ
765 : : end do ! idat
766 : :
767 : : #else
768 [ + - ]: 191868 : if(gpu_option==ABI_GPU_DISABLED) then
769 : : !$OMP PARALLEL DO COLLAPSE(2)
770 [ + + ]: 383736 : do idat=1,ndat
771 [ + + ]: 1244635 : do idat_occ=1,ndat_occ
772 [ + + ]: 3870429685 : do ifft=1,nfftf
773 : 3869376918 : rhog_munu(1,ifft,idat_occ,idat) = rhog_munu(1,ifft,idat_occ,idat) * vqg(ifft)
774 : 3870237817 : rhog_munu(2,ifft,idat_occ,idat) = rhog_munu(2,ifft,idat_occ,idat) * vqg(ifft)
775 : : end do
776 : : end do ! idat_occ
777 : : end do ! idat
778 : : else if(gpu_option==ABI_GPU_OPENMP) then
779 : : #ifdef HAVE_OPENMP_OFFLOAD
780 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(2) &
781 : : !$OMP& MAP(to:rhog_munu,vqg) PRIVATE(idat,idat_occ)
782 : : do idat=1,ndat
783 : : do idat_occ=1,ndat_occ
784 : : !$OMP PARALLEL DO PRIVATE(ifft)
785 : : do ifft=1,nfftf
786 : : rhog_munu(1,ifft,idat_occ,idat) = rhog_munu(1,ifft,idat_occ,idat) * vqg(ifft)
787 : : rhog_munu(2,ifft,idat_occ,idat) = rhog_munu(2,ifft,idat_occ,idat) * vqg(ifft)
788 : : end do
789 : : end do ! idat_occ
790 : : end do ! idat
791 : : #endif
792 : : end if ! gpu_option
793 : :
794 : 191868 : call timab(1515,2,tsec) ; call timab(1513,-1,tsec) ; call timab(1545,-2,tsec)
795 : 191868 : call fourdp(cplex_fock,rhog_munu,vfock,+1,mpi_enreg,nfftf,ndat*ndat_occ,ngfftf,tim_fourdp_fock_getghc,gpu_option=gpu_option)
796 : 191868 : call timab(1513,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1545,-1,tsec)
797 : : #endif
798 : 191868 : call timab(1525,-2,tsec) ; call timab(1545,-2,tsec)
799 : :
800 : : !===============================================================
801 : : !======== Calculate Dij_Fock_hat contribution in case of PAW ===
802 : : !===============================================================
803 : :
804 : 191868 : call timab(1526,-1,tsec) ; call timab(1546,-1,tsec)
805 : :
806 [ + + ]: 191868 : if (fockcommon%usepaw==1) then
807 [ + + ]: 40440 : qphon=qvec_j;nfftotf=product(ngfftf(1:3))
808 : 10110 : ndij=nspden_fock
809 : : ! dimekb1 is dimensioned as cplex_dij*lmnmax*(lmnmax+1)/2
810 : 10110 : cplex_dij=2*gs_ham%dimekb1/(gs_ham%lmnmax*(gs_ham%lmnmax+1))
811 [ - + + - : 70770 : ABI_MALLOC(dijhat,(gs_ham%dimekb1,natom,ndij,ndat_tot,cplex_fock))
- + + - -
+ + - - +
+ - - + +
- + - + -
- + - + -
+ ]
812 [ + + + + : 6125682 : dijhat=zero
+ + + + +
+ ]
813 : :
814 : : #ifdef HAVE_OPENMP_OFFLOAD
815 : : !!$OMP TARGET UPDATE FROM(vfock) IF(gpu_option==ABI_GPU_OPENMP)
816 : : #endif
817 : 10110 : iatm=0
818 [ + + ]: 20934 : do itypat=1,gs_ham%ntypat
819 : 10824 : lmn2_size=fockcommon%pawtab(itypat)%lmn2_size
820 [ - + + - : 54120 : ABI_MALLOC(dijhat_tmp,(cplex_fock*cplex_dij*lmn2_size,ndij*ndat_occ*ndat,gs_ham%nattyp(itypat)))
- + + - -
+ + - + -
+ - - + -
+ - + ]
821 : : call pawdijhat_ndat(dijhat_tmp,cplex_dij,cplex_fock,gs_ham%gprimd,iatm,&
822 : : & natom,ndij,nfftf,nfftotf,nspden_fock,nspden_fock,ndat_occ*ndat,&
823 : : & gs_ham%nattyp(itypat),fockbz%pawang,fockcommon%pawfgrtab,&
824 : : & fockcommon%pawtab(itypat),vfock,qphon,gs_ham%ucvol,gs_ham%xred,&
825 [ + - - - : 10824 : & gpu_option=gpu_option)
- - - - -
+ ]
826 : : !$OMP PARALLEL DO COLLAPSE(2) PRIVATE(ia,idat,idat_occ,idat_tot,ii,iatom,ind)
827 [ + + ]: 31044 : do ia=1,gs_ham%nattyp(itypat)
828 [ + + ]: 51264 : do idat=1,ndat
829 [ + + ]: 120642 : do idat_occ=1,ndat_occ
830 : 80202 : idat_tot=(idat-1)*ndat_occ+idat_occ
831 : 80202 : iatom=iatm+ia
832 [ + + ]: 260826 : do ii=1,cplex_fock
833 : 160404 : ind=(ii-1)*lmn2_size*cplex_dij
834 : : dijhat(1:cplex_dij*lmn2_size,iatom,:,idat_tot,ii)=&
835 : : & dijhat_tmp(ind+1:ind+cplex_dij*lmn2_size,&
836 [ + + + + ]: 6027912 : 1+(idat_occ-1)*ndij+(idat-1)*ndat_occ*ndij:idat_occ*ndij+(idat-1)*ndat_occ*ndij,ia)
837 : : end do
838 : : end do ! idat_occ
839 : : end do ! idat
840 : : end do ! ia
841 [ - + ]: 10824 : ABI_FREE(dijhat_tmp)
842 : 20934 : iatm=iatm+gs_ham%nattyp(itypat)
843 : : end do
844 : :
845 : : !*NOTE Build a replica of cwaveocc_prj covering the full ndat*ndat_occ batch: cwaveocc_prj
846 : : !*itself does not depend on idat, but gemm_nonlop requires cprjin to be sized on the
847 : : !*full batch passed as its "ndat" argument.
848 : : !*This is a bit tedious but remains cheaper than calling nonlop ndat times
849 [ - + + - : 160743 : ABI_MALLOC(cwaveocc_prj_rep,(natom,nspinor*ndat_tot))
- + + - +
- + - - +
- + + + +
+ ]
850 [ + - + + ]: 30330 : call pawcprj_alloc(cwaveocc_prj_rep,cwaveocc_prj(1,1)%ncpgr,cwaveocc_prj(:,1)%nlmn)
851 [ + + ]: 20220 : do idat=1,ndat
852 : : call pawcprj_copy(cwaveocc_prj,&
853 : 20220 : & cwaveocc_prj_rep(:,(idat-1)*nspinor*ndat_occ+1:idat*nspinor*ndat_occ))
854 : : end do
855 : :
856 : 10110 : signs=2; cpopt=2;idir=0; paw_opt=1;nnlout=1;tim_nonlop=17
857 : :
858 [ + + ]: 10110 : if(need_ghc) then
859 : 7952 : choice=1
860 : 7952 : call timab(1515,2,tsec) ; call timab(1514,-1,tsec) ; call timab(1546,-2,tsec)
861 : : call nonlop(choice,cpopt,cwaveocc_prj_rep,enlout_dum,gs_ham,idir,(/zero/),&
862 : : & mpi_enreg,ndat_tot,nnlout,paw_opt,signs,gsc_dum,tim_nonlop,vectin_dum,&
863 : : & gvnlxc,enl_ndat=dijhat,&
864 : 7952 : & select_k=K_H_KPRIME)
865 : :
866 [ + - ]: 7952 : if(gpu_option==ABI_GPU_DISABLED) then
867 [ + + ]: 15904 : do idat=1,ndat
868 [ + + ]: 47368 : do idat_occ=1,ndat_occ
869 : 31464 : idat_tot=(idat-1)*ndat_occ+idat_occ
870 : : ghc2(:,1+(idat-1)*npw*nspinor:idat*npw*nspinor) = ghc2(:,1+(idat-1)*npw*nspinor:idat*npw*nspinor)&
871 [ + + + + ]: 24516644 : & - gvnlxc(:,1+(idat_tot-1)*npw*nspinor:idat_tot*npw*nspinor) * occ(idat_occ) * wtk
872 : : end do ! idat_occ
873 : : end do ! idat
874 : : else if(gpu_option==ABI_GPU_OPENMP) then
875 : : #ifdef HAVE_OPENMP_OFFLOAD
876 : : !$OMP TARGET TEAMS DISTRIBUTE PARALLEL DO COLLAPSE(2) &
877 : : !$OMP& MAP(to:gvnlxc,ghc2,occ) PRIVATE(idat,ipw,idat_occ,idat_tot)
878 : : do idat=1,ndat
879 : : do ipw=1,npw
880 : : do idat_occ=1,ndat_occ
881 : : idat_tot=(idat-1)*ndat_occ+idat_occ
882 : : ghc2(1,ipw+(idat-1)*npw*nspinor) = ghc2(1,ipw+(idat-1)*npw*nspinor)&
883 : : & - gvnlxc(1,ipw+(idat_tot-1)*npw*nspinor) * occ(idat_occ) * wtk
884 : : ghc2(2,ipw+(idat-1)*npw*nspinor) = ghc2(2,ipw+(idat-1)*npw*nspinor)&
885 : : & - gvnlxc(2,ipw+(idat_tot-1)*npw*nspinor) * occ(idat_occ) * wtk
886 : : end do
887 : : end do
888 : : end do ! idat
889 : : #endif
890 : : end if
891 : 7952 : call timab(1514,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1546,-1,tsec)
892 : : end if
893 : :
894 : : ! Forces calculation
895 : :
896 [ + + + - ]: 10110 : if (fockcommon%optfor.and.(fockcommon%ieigen/=0)) then
897 [ - + + - : 17976 : ABI_MALLOC(vdotr,(ndat_occ,ndat,3,natom))
- + + - +
- - + + -
+ - + - -
+ - + ]
898 [ - + + - : 8988 : ABI_MALLOC(vdoti,(ndat_tot))
+ - + - -
+ - + ]
899 [ - + + - : 14980 : ABI_MALLOC(for1,(ndat_occ,3,natom,ndat))
- + + - -
+ + - + -
+ - - + -
+ ]
900 : : #ifdef HAVE_OPENMP_OFFLOAD
901 : : !$OMP TARGET ENTER DATA MAP(alloc:vdotr,vdoti,for1) IF(gpu_option==ABI_GPU_OPENMP)
902 : : #endif
903 : : !*NOTE:cwavef_rep replicates cwavef ndat_occ times (per idat block) so that single
904 : : !*nonlop and dotprod_g_batch_full call can process the whole ndat*ndat_occ batch at once.
905 [ - + + - : 8988 : ABI_MALLOC(cwavef_rep,(2,npw*nspinor*ndat_tot))
+ - + - -
+ - + ]
906 : : #ifdef HAVE_OPENMP_OFFLOAD
907 : : !$OMP TARGET ENTER DATA MAP(alloc:cwavef_rep) IF(gpu_option==ABI_GPU_OPENMP)
908 : : #endif
909 [ + - ]: 2996 : if(gpu_option==ABI_GPU_DISABLED) then
910 [ + + ]: 5992 : do idat=1,ndat
911 [ + + ]: 17632 : do idat_occ=1,ndat_occ
912 : 11640 : idat_tot=(idat-1)*ndat_occ+idat_occ
913 : : cwavef_rep(:,1+(idat_tot-1)*npw*nspinor:idat_tot*npw*nspinor)=&
914 [ + + + + ]: 17584004 : & cwavef(:,1+(idat-1)*npw*nspinor:idat*npw*nspinor)
915 : : end do ! idat_occ
916 : : end do ! idat
917 : : else if(gpu_option==ABI_GPU_OPENMP) then
918 : : #ifdef HAVE_OPENMP_OFFLOAD
919 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(2) MAP(to:cwavef_rep,cwavef) PRIVATE(idat,idat_occ,idat_tot)
920 : : do idat=1,ndat
921 : : do idat_occ=1,ndat_occ
922 : : idat_tot=(idat-1)*ndat_occ+idat_occ
923 : : !$OMP PARALLEL DO PRIVATE(ipw)
924 : : do ipw=1,npw*nspinor
925 : : cwavef_rep(1,ipw+(idat_tot-1)*npw*nspinor)=cwavef(1,ipw+(idat-1)*npw*nspinor)
926 : : cwavef_rep(2,ipw+(idat_tot-1)*npw*nspinor)=cwavef(2,ipw+(idat-1)*npw*nspinor)
927 : : end do
928 : : end do ! idat_occ
929 : : end do ! idat
930 : : #endif
931 : : end if
932 [ + + + + : 114780 : choice=2; vdotr=zero;doti=zero;cpopt=4;tim_nonlop=17
+ + + + ]
933 [ + + ]: 11984 : do idir=1,3
934 [ + + ]: 29960 : do iatom=1,natom
935 : 17976 : call timab(1515,2,tsec) ; call timab(1514,-1,tsec) ; call timab(1546,-2,tsec)
936 : : call nonlop(choice,cpopt,cwaveocc_prj_rep,enlout_dum,gs_ham,idir,(/zero/),mpi_enreg,&
937 : : & ndat_tot,nnlout,paw_opt,signs,gsc_dum,tim_nonlop,vectin_dum,&
938 : : & forout,enl_ndat=dijhat,iatom_only=iatom,&
939 : 17976 : & select_k=K_H_KPRIME)
940 : 17976 : call timab(1514,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1546,-1,tsec)
941 : : call dotprod_g_batch_full(vdotr(:,:,idir,iatom),vdoti,gs_ham%istwf_k,npw,ndat_tot,2,&
942 : 26964 : & cwavef_rep,forout,mpi_enreg%me_g0,mpi_enreg%comm_fft,gpu_option=gpu_option)
943 : : end do ! iatom
944 : : end do ! idir
945 : : #ifdef HAVE_OPENMP_OFFLOAD
946 : : !$OMP TARGET EXIT DATA MAP(delete:cwavef_rep) IF(gpu_option==ABI_GPU_OPENMP)
947 : : #endif
948 [ - + ]: 2996 : ABI_FREE(cwavef_rep)
949 : :
950 [ + - ]: 2996 : if(gpu_option==ABI_GPU_DISABLED) then
951 [ + + ]: 5992 : do idat=1,ndat
952 [ + + ]: 11984 : do iatom=1,natom
953 [ + + ]: 26964 : do idir=1,3
954 [ + + ]: 93808 : do idat_occ=1,ndat_occ
955 : 69840 : esum=0
956 [ + + ]: 59108256 : do ifft=1,atom_nfgd(iatom)
957 : 59038416 : ind=atom_ifftsph(ifft,iatom)
958 : : esum = esum &
959 : : & + vfock(2*ind-1,idat_occ,idat)*grnhat_12(1,ifft,1,idir,iatom,idat_occ,idat) &
960 : 59108256 : & - vfock(2*ind,idat_occ,idat)*grnhat_12(2,ifft,1,idir,iatom,idat_occ,idat)
961 : : end do
962 : 87816 : for1(idat_occ,idir,iatom,idat)=esum
963 : : end do ! idat_occ
964 : : end do ! idir
965 : : end do ! iatom
966 : : end do ! idat
967 : : else if(gpu_option==ABI_GPU_OPENMP) then
968 : : #ifdef HAVE_OPENMP_OFFLOAD
969 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(2) &
970 : : !$OMP& MAP(to:vfock,grnhat_12,for1,atom_nfgd,atom_ifftsph) &
971 : : !$OMP& PRIVATE(ifft,ind,iatom) PRIVATE(idat_occ,idir,esum)
972 : : do idat=1,ndat
973 : : do iatom=1,natom
974 : : !$OMP PARALLEL DO COLLAPSE(2) PRIVATE(esum,ind,idir,idat_occ,ifft)
975 : : do idir=1,3
976 : : do idat_occ=1,ndat_occ
977 : : esum=0
978 : : do ifft=1,atom_nfgd(iatom)
979 : : ind=atom_ifftsph(ifft,iatom)
980 : : esum = esum &
981 : : & + vfock(2*ind-1,idat_occ,idat)*grnhat_12(1,ifft,1,idir,iatom,idat_occ,idat) &
982 : : & - vfock(2*ind,idat_occ,idat)*grnhat_12(2,ifft,1,idir,iatom,idat_occ,idat)
983 : : end do
984 : : for1(idat_occ,idir,iatom,idat)=esum
985 : : end do ! idat_occ
986 : : end do ! idir
987 : : end do ! iatom
988 : : end do ! idat
989 : : !$OMP TARGET UPDATE FROM(for1)
990 : : #endif
991 : : end if
992 : :
993 : : if(.true.) then
994 : : !if(gpu_option==ABI_GPU_DISABLED) then
995 : : #ifdef HAVE_OPENMP_OFFLOAD
996 : : !$OMP TARGET UPDATE FROM(vdotr) IF(gpu_option==ABI_GPU_OPENMP)
997 : : #endif
998 [ + + ]: 5992 : do idat=1,ndat
999 [ + + ]: 11984 : do iatom=1,natom
1000 [ + + ]: 26964 : do idir=1,3
1001 [ + + ]: 93808 : do idat_occ=1,ndat_occ
1002 : : for12(idir)=rprimd(1,idir)*for1(idat_occ,1,iatom,idat)&
1003 : : & +rprimd(2,idir)*for1(idat_occ,2,iatom,idat)&
1004 : 69840 : & +rprimd(3,idir)*for1(idat_occ,3,iatom,idat)
1005 : : forikpt(idir,iatom,idat)=forikpt(idir,iatom,idat)&
1006 : : & - (for12(idir)*gs_ham%ucvol/nfftf+vdotr(idat_occ,idat,idir,iatom))&
1007 : 87816 : & * occ(idat_occ) * wtk
1008 : : end do ! idat_occ
1009 : : end do ! idir
1010 : : end do ! iatom
1011 : : end do ! idat
1012 : : end if
1013 : : #ifdef HAVE_OPENMP_OFFLOAD
1014 : : !$OMP TARGET EXIT DATA MAP(delete:vdotr,vdoti,for1) IF(gpu_option==ABI_GPU_OPENMP)
1015 : : #endif
1016 [ - + ]: 2996 : ABI_FREE(vdotr)
1017 [ - + ]: 2996 : ABI_FREE(vdoti)
1018 [ - + ]: 2996 : ABI_FREE(for1)
1019 : : end if
1020 : :
1021 : : ! Stresses calculation
1022 [ + + + - ]: 10110 : if (fockcommon%optstr.and.(fockcommon%ieigen/=0)) then
1023 : 1870 : signs=2;choice=3;cpopt=4;tim_nonlop=17
1024 : :
1025 : : ! first contribution
1026 [ - + + - : 9350 : ABI_MALLOC(vdotr,(ndat_occ,ndat,6,1))
- + + - +
- + - + -
- + - + ]
1027 [ - + + - : 5610 : ABI_MALLOC(vdoti,(ndat_tot))
+ - + - -
+ - + ]
1028 : : #ifdef HAVE_OPENMP_OFFLOAD
1029 : : !$OMP TARGET ENTER DATA MAP(alloc:vdotr,vdoti) IF(gpu_option==ABI_GPU_OPENMP)
1030 : : #endif
1031 : : !*cwavef_rep replicates cwavef ndat_occ times (per idat block) so that a single
1032 : : !*dotprod_g_batch_full call can process the whole ndat*ndat_occ batch at once,
1033 : : !*matching the merged nonlop call below.
1034 [ - + + - : 5610 : ABI_MALLOC(cwavef_rep,(2,npw*nspinor*ndat_tot))
+ - + - -
+ - + ]
1035 : : #ifdef HAVE_OPENMP_OFFLOAD
1036 : : !$OMP TARGET ENTER DATA MAP(alloc:cwavef_rep) IF(gpu_option==ABI_GPU_OPENMP)
1037 : : #endif
1038 [ + - ]: 1870 : if(gpu_option==ABI_GPU_DISABLED) then
1039 [ + + ]: 3740 : do idat=1,ndat
1040 [ + + ]: 11225 : do idat_occ=1,ndat_occ
1041 : 7485 : idat_tot=(idat-1)*ndat_occ+idat_occ
1042 : : cwavef_rep(:,1+(idat_tot-1)*npw*nspinor:idat_tot*npw*nspinor)=&
1043 [ + + + + ]: 4097380 : & cwavef(:,1+(idat-1)*npw*nspinor:idat*npw*nspinor)
1044 : : end do ! idat_occ
1045 : : end do ! idat
1046 : : else if(gpu_option==ABI_GPU_OPENMP) then
1047 : : #ifdef HAVE_OPENMP_OFFLOAD
1048 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(2) MAP(to:cwavef_rep,cwavef) PRIVATE(idat,idat_occ,idat_tot)
1049 : : do idat=1,ndat
1050 : : do idat_occ=1,ndat_occ
1051 : : idat_tot=(idat-1)*ndat_occ+idat_occ
1052 : : !$OMP PARALLEL DO PRIVATE(ipw)
1053 : : do ipw=1,npw*nspinor
1054 : : cwavef_rep(1,ipw+(idat_tot-1)*npw*nspinor)=cwavef(1,ipw+(idat-1)*npw*nspinor)
1055 : : cwavef_rep(2,ipw+(idat_tot-1)*npw*nspinor)=cwavef(2,ipw+(idat-1)*npw*nspinor)
1056 : : end do
1057 : : end do ! idat_occ
1058 : : end do ! idat
1059 : : #endif
1060 : : end if
1061 [ + + ]: 13090 : do idir=1,6
1062 : 11220 : call timab(1515,2,tsec) ; call timab(1514,-1,tsec) ; call timab(1546,-2,tsec)
1063 : : call nonlop(choice,cpopt,cwaveocc_prj_rep,enlout_dum,gs_ham,idir,(/zero/),mpi_enreg,&
1064 : : & ndat_tot,nnlout,paw_opt,signs,gsc_dum,tim_nonlop,vectin_dum,&
1065 : 11220 : & strout,enl_ndat=dijhat,select_k=K_H_KPRIME)
1066 : 11220 : call timab(1514,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1546,-1,tsec)
1067 : : call dotprod_g_batch_full(vdotr(:,:,idir,1),vdoti,gs_ham%istwf_k,npw,ndat_tot,2,&
1068 : : & cwavef_rep,strout,&
1069 : 11220 : & mpi_enreg%me_g0,mpi_enreg%comm_fft,gpu_option=gpu_option)
1070 : : #ifdef HAVE_OPENMP_OFFLOAD
1071 : : !$OMP TARGET UPDATE FROM(vdotr) IF(gpu_option==ABI_GPU_OPENMP)
1072 : : #endif
1073 [ + + ]: 24310 : do idat=1,ndat
1074 [ + + ]: 67350 : do idat_occ=1,ndat_occ
1075 : : fockcommon%stress_ikpt(idir,fockcommon%ieigen+idat-1)=fockcommon%stress_ikpt(idir,fockcommon%ieigen+idat-1)-&
1076 : 56130 : & vdotr(idat_occ,idat,idir,1)*occ(idat_occ)*wtk/gs_ham%ucvol
1077 : : end do ! idat_occ
1078 : : end do ! idat
1079 : : end do ! idir
1080 : : #ifdef HAVE_OPENMP_OFFLOAD
1081 : : !$OMP TARGET EXIT DATA MAP(delete:vdotr,vdoti,cwavef_rep) IF(gpu_option==ABI_GPU_OPENMP)
1082 : : #endif
1083 [ - + ]: 1870 : ABI_FREE(vdotr)
1084 [ - + ]: 1870 : ABI_FREE(vdoti)
1085 [ - + ]: 1870 : ABI_FREE(cwavef_rep)
1086 : :
1087 : 1870 : stress_ikpt => fockcommon%stress_ikpt
1088 : 1870 : ieigen = fockcommon%ieigen
1089 : :
1090 : : ! second contribution
1091 : : !if(.true.) then
1092 [ + - ]: 1870 : if(gpu_option==ABI_GPU_DISABLED) then
1093 : : #ifdef HAVE_OPENMP_OFFLOAD
1094 : : !$OMP TARGET UPDATE FROM(vfock,grnhat_12) IF(gpu_option==ABI_GPU_OPENMP)
1095 : : #endif
1096 [ - + + - : 7480 : ABI_MALLOC(strdat, (3,3,ndat_occ,ndat))
- + + - +
- + - - +
- + ]
1097 [ + + + + : 101045 : strdat=zero
+ + + + ]
1098 [ + + ]: 3740 : do idat=1,ndat
1099 [ + + ]: 11225 : do idat_occ=1,ndat_occ
1100 [ + + ]: 31810 : do idir=1,3
1101 [ + + ]: 97305 : do idir1=1,3
1102 : 67365 : esum=0
1103 [ + + ]: 202095 : do iatom=1,natom
1104 [ + + ]: 104707395 : do ifft=1,atom_nfgd(iatom)
1105 : 104505300 : ind = atom_ifftsph(ifft,iatom)
1106 : : esum = esum + (vfock(2*ind-1,idat_occ,idat)*grnhat_12(1,ifft,1,idir,iatom,idat_occ,idat)&
1107 : : & - vfock(2*ind,idat_occ,idat)*grnhat_12(2,ifft,1,idir,iatom,idat_occ,idat))&
1108 : 104640030 : & * atom_rfgd(idir1,ifft,iatom)
1109 : : end do
1110 : : end do
1111 : 89820 : strdat(idir,idir1,idat_occ,idat)=esum
1112 : : end do
1113 : : end do
1114 : : end do ! idat_occ
1115 : : end do ! idat
1116 [ + + ]: 3740 : do idat=1,ndat
1117 [ + + ]: 11225 : do idat_occ=1,ndat_occ
1118 [ + + ]: 29940 : do idir=1,3
1119 : 29940 : fockstr(idir)=strdat(idir,idir,idat_occ,idat)
1120 : : end do
1121 : 7485 : fockstr(4)=(strdat(3,2,idat_occ,idat)+strdat(2,3,idat_occ,idat))*half
1122 : 7485 : fockstr(5)=(strdat(3,1,idat_occ,idat)+strdat(1,3,idat_occ,idat))*half
1123 : 7485 : fockstr(6)=(strdat(1,2,idat_occ,idat)+strdat(2,1,idat_occ,idat))*half
1124 [ + + ]: 54265 : do idir=1,6
1125 : : fockcommon%stress_ikpt(idir,fockcommon%ieigen+idat-1)=fockcommon%stress_ikpt(idir,fockcommon%ieigen+idat-1)+&
1126 : 52395 : & fockstr(idir)/nfftf*occ(idat_occ)*wtk
1127 : : end do
1128 : : end do ! idat_occ
1129 : : end do ! idat
1130 [ - + ]: 1870 : ABI_FREE(strdat)
1131 : : else if(gpu_option==ABI_GPU_OPENMP) then
1132 : : #ifdef HAVE_OPENMP_OFFLOAD
1133 : : ABI_MALLOC(strdat, (3,3,ndat_occ,ndat))
1134 : : strdat=zero
1135 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(4) MAP(tofrom:strdat) &
1136 : : !$OMP& MAP(to:vfock,grnhat_12,atom_nfgd,atom_rfgd,atom_ifftsph) &
1137 : : !$OMP& PRIVATE(idat,idat_occ,idir,idir1,esum,iatom)
1138 : : do idat=1,ndat
1139 : : do idat_occ=1,ndat_occ
1140 : : do idir=1,3
1141 : : do idir1=1,3
1142 : : do iatom=1,natom
1143 : : esum=0
1144 : : !$OMP PARALLEL DO PRIVATE(ifft,ind) REDUCTION(+:esum)
1145 : : do ifft=1,atom_nfgd(iatom)
1146 : : ind=atom_ifftsph(ifft,iatom)
1147 : : esum = esum + (vfock(2*ind-1,idat_occ,idat)*grnhat_12(1,ifft,1,idir,iatom,idat_occ,idat)&
1148 : : & - vfock(2*ind,idat_occ,idat)*grnhat_12(2,ifft,1,idir,iatom,idat_occ,idat))&
1149 : : & * atom_rfgd(idir1,ifft,iatom)
1150 : : end do
1151 : : strdat(idir,idir1,idat_occ,idat)=strdat(idir,idir1,idat_occ,idat)+esum
1152 : : end do
1153 : : end do
1154 : : end do
1155 : : end do ! idat_occ
1156 : : end do ! idat
1157 : : do idat=1,ndat
1158 : : do idat_occ=1,ndat_occ
1159 : : do idir=1,3
1160 : : fockstr(idir)=strdat(idir,idir,idat_occ,idat)
1161 : : end do
1162 : : fockstr(4)=(strdat(3,2,idat_occ,idat)+strdat(2,3,idat_occ,idat))*half
1163 : : fockstr(5)=(strdat(3,1,idat_occ,idat)+strdat(1,3,idat_occ,idat))*half
1164 : : fockstr(6)=(strdat(1,2,idat_occ,idat)+strdat(2,1,idat_occ,idat))*half
1165 : : do idir=1,6
1166 : : stress_ikpt(idir,ieigen+idat-1)=stress_ikpt(idir,ieigen+idat-1)+&
1167 : : & fockstr(idir)/nfftf*occ(idat_occ)*wtk
1168 : : end do
1169 : : end do ! idat_occ
1170 : : end do ! idat
1171 : : ABI_FREE(strdat)
1172 : : #endif
1173 : : end if
1174 : :
1175 : : ! third contribution
1176 [ + - ]: 1870 : if(gpu_option==ABI_GPU_DISABLED) then
1177 [ + + ]: 3740 : do idat=1,ndat
1178 [ + + ]: 11225 : do idat_occ=1,ndat_occ
1179 : 7485 : doti=zero
1180 [ + + ]: 22455 : do iatom=1,natom
1181 [ + + ]: 11634155 : do ifft=1,atom_nfgd(iatom)
1182 : 11611700 : ind = atom_ifftsph(ifft,iatom)
1183 : : doti = doti &
1184 : : & + vfock(2*ind-1,idat_occ,idat) * rho12(1,ifft,nspinor,idat_occ,idat,iatom)&
1185 : 11626670 : & - vfock(2*ind,idat_occ,idat) * rho12(2,ifft,nspinor,idat_occ,idat,iatom)
1186 : : end do
1187 : : end do
1188 [ + + ]: 31810 : fockcommon%stress_ikpt(1:3,fockcommon%ieigen+idat-1)=fockcommon%stress_ikpt(1:3,fockcommon%ieigen+idat-1)-doti/nfftf*occ(idat_occ)*wtk
1189 : : end do ! idat_occ
1190 : : end do ! idat
1191 : : else if(gpu_option==ABI_GPU_OPENMP) then
1192 : : #ifdef HAVE_OPENMP_OFFLOAD
1193 : : ABI_MALLOC(strdat, (1,1,ndat_occ,ndat))
1194 : : strdat=zero
1195 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(2) &
1196 : : !$OMP& MAP(tofrom:strdat) MAP(to:vfock,rho12,atom_nfgd,atom_ifftsph) &
1197 : : !$OMP& PRIVATE(idat,idat_occ,idir,idir1,doti)
1198 : : do idat=1,ndat
1199 : : do idat_occ=1,ndat_occ
1200 : : doti=zero
1201 : : !$OMP PARALLEL DO PRIVATE(iatom,ifft,ind) REDUCTION(+:doti)
1202 : : do iatom=1,natom
1203 : : do ifft=1,atom_nfgd(iatom)
1204 : : ind = atom_ifftsph(ifft,iatom)
1205 : : doti = doti &
1206 : : & + vfock(2*ind-1,idat_occ,idat) * rho12(1,ifft,nspinor,idat_occ,idat,iatom)&
1207 : : & - vfock(2*ind,idat_occ,idat) * rho12(2,ifft,nspinor,idat_occ,idat,iatom)
1208 : : end do
1209 : : end do
1210 : : strdat(1,1,idat_occ,idat)=doti
1211 : : end do ! idat_occ
1212 : : end do ! idat
1213 : : do idat=1,ndat
1214 : : do idat_occ=1,ndat_occ
1215 : : fockcommon%stress_ikpt(1:3,fockcommon%ieigen+idat-1)=fockcommon%stress_ikpt(1:3,fockcommon%ieigen+idat-1) &
1216 : : -strdat(1,1,idat_occ,idat)/nfftf*occ(idat_occ)*wtk
1217 : : end do ! idat_occ
1218 : : end do ! idat
1219 : : ABI_FREE(strdat)
1220 : : #endif
1221 : : end if
1222 : : end if ! end stresses
1223 : :
1224 : 10110 : call pawcprj_free(cwaveocc_prj_rep)
1225 [ - + ]: 10110 : ABI_FREE(cwaveocc_prj_rep)
1226 [ - + ]: 10110 : ABI_FREE(dijhat)
1227 : : end if !end PAW
1228 : 191868 : call timab(1526,2,tsec) ; call timab(1546,-2,tsec)
1229 : :
1230 : : ! =============================================================
1231 : : ! === Apply the local potential vfockloc_munu to cwaveocc_r ===
1232 : : ! =============================================================
1233 : 191868 : call timab(1527,-1,tsec)
1234 [ + - ]: 191868 : if(gpu_option==ABI_GPU_DISABLED) then
1235 : : ABI_UNUSED(esumi)
1236 : : !$OMP PARALLEL DO &
1237 : : !$OMP& PRIVATE(ind,recwocc,imcwocc,revloc,imvloc)
1238 [ + + ]: 383736 : do idat=1,ndat
1239 [ + + ]: 1244635 : do idat_occ=1,ndat_occ
1240 [ + + ]: 14678740 : do i3=1,ngfftf(3)
1241 [ + + ]: 236443177 : do i2=1,ngfftf(2)
1242 [ + + ]: 4104959196 : do i1=1,ngfftf(1)
1243 : 3869376918 : ind=i1+(i2-1)*ngfftf(1)+(i3-1)*ngfftf(2)*ngfftf(1)
1244 : 3869376918 : revloc=vfock(2*ind-1,idat_occ,idat) ; imvloc=vfock(2*ind,idat_occ,idat)
1245 : 3869376918 : recwocc=cwaveocc_r(1,i1,i2,i3,idat_occ)
1246 : 3869376918 : imcwocc=cwaveocc_r(2,i1,i2,i3,idat_occ)
1247 : 3869376918 : vlocpsi_r(2*ind-1,idat)=vlocpsi_r(2*ind-1,idat)-(revloc*recwocc-imvloc*imcwocc)*occ(idat_occ)*wtk
1248 : 4091333223 : vlocpsi_r(2*ind ,idat)=vlocpsi_r(2*ind ,idat)-(revloc*imcwocc+imvloc*recwocc)*occ(idat_occ)*wtk
1249 : : end do
1250 : : end do
1251 : : end do
1252 : : end do ! idat_occ
1253 : : end do ! idat
1254 : : else if(gpu_option==ABI_GPU_OPENMP) then
1255 : : #ifdef HAVE_OPENMP_OFFLOAD
1256 : : !$OMP TARGET TEAMS DISTRIBUTE COLLAPSE(3) &
1257 : : !$OMP& MAP(to:vlocpsi_r,cwaveocc_r,occ,vfock) PRIVATE(idat,i3,i2)
1258 : : do idat=1,ndat
1259 : : do i3=1,ngfftf(3)
1260 : : do i2=1,ngfftf(2)
1261 : : !$OMP PARALLEL DO PRIVATE(ind,i1,idat_occ,recwocc,imcwocc,revloc,imvloc,esum,esumi)
1262 : : do i1=1,ngfftf(1)
1263 : : ind=i1+(i2-1)*ngfftf(1)+(i3-1)*ngfftf(2)*ngfftf(1)
1264 : : esum=zero ; esumi=zero
1265 : : do idat_occ=1,ndat_occ
1266 : : revloc=vfock(2*ind-1,idat_occ,idat) ; imvloc=vfock(2*ind,idat_occ,idat)
1267 : : recwocc=cwaveocc_r(1,i1,i2,i3,idat_occ)
1268 : : imcwocc=cwaveocc_r(2,i1,i2,i3,idat_occ)
1269 : : esum=esum+(revloc*recwocc-imvloc*imcwocc)*occ(idat_occ)
1270 : : esumi=esumi+(revloc*imcwocc+imvloc*recwocc)*occ(idat_occ)
1271 : : end do
1272 : : vlocpsi_r(2*ind-1,idat)=vlocpsi_r(2*ind-1,idat)-esum*wtk
1273 : : vlocpsi_r(2*ind ,idat)=vlocpsi_r(2*ind ,idat)-esumi*wtk
1274 : : end do
1275 : : end do
1276 : : end do
1277 : : end do ! idat
1278 : : #endif
1279 : : end if
1280 [ + - ]: 191868 : if (allocated(fockbz%cgocc)) then
1281 : : #ifdef HAVE_OPENMP_OFFLOAD
1282 : : !$OMP TARGET EXIT DATA MAP(delete:cwaveocc_r) IF(gpu_option==ABI_GPU_OPENMP)
1283 : : #endif
1284 [ - + ]: 191868 : ABI_FREE(cwaveocc_r)
1285 : : end if
1286 : 394826 : call timab(1527,2,tsec)
1287 : : end do ! jband
1288 : :
1289 : : ! ========================================================
1290 : : ! === End of loop : update of shifts and deallocations ===
1291 : : ! ==============================:=========================
1292 : : !* Update of the shifts to be applied (reminder : mkmem is not 0, nspinor=1)
1293 : 202958 : call timab(1528,1,tsec)
1294 : 202958 : jcg=jcg+npwj*nband_k
1295 : 202958 : jbg=jbg+nband_k
1296 : 202958 : bdtot_jindex=bdtot_jindex+nband_k
1297 : 202958 : jkg=jkg+npwj
1298 [ + + ]: 202958 : if (fockcommon%usepaw==1) then
1299 [ - + ]: 10110 : ABI_FREE(enlout_dum)
1300 [ - + ]: 10110 : ABI_FREE(vectin_dum)
1301 [ - + ]: 10110 : ABI_FREE(ffnl_kp_dum)
1302 [ - + ]: 10110 : ABI_FREE(kpg_kp)
1303 : : end if
1304 [ + - ]: 202958 : if (associated(gs_ham%ph3d_kp)) then
1305 : 202958 : ABI_FREE(gs_ham%ph3d_kp)
1306 : : end if
1307 : : #ifdef HAVE_OPENMP_OFFLOAD
1308 : : !$OMP TARGET EXIT DATA MAP(delete:rhor_munu,rhog_munu,vfock,occ,rho12,gvnlxc) IF(gpu_option==ABI_GPU_OPENMP)
1309 : : #endif
1310 [ - + ]: 202958 : ABI_FREE(rhor_munu)
1311 [ - + ]: 202958 : ABI_FREE(rhog_munu)
1312 [ - + ]: 202958 : ABI_FREE(vfock)
1313 [ - + ]: 202958 : ABI_FREE(occ)
1314 : : !*Additional arrays in case of paw
1315 [ + + ]: 202958 : if (fockcommon%usepaw==1) then
1316 [ + + + - ]: 10110 : if ((fockcommon%optfor).and.(fockcommon%ieigen/=0)) then
1317 : : #ifdef HAVE_OPENMP_OFFLOAD
1318 : : !$OMP TARGET EXIT DATA MAP(delete:forout) IF(gpu_option==ABI_GPU_OPENMP)
1319 : : #endif
1320 [ - + ]: 2996 : ABI_FREE(forout)
1321 : : end if
1322 [ + + + - ]: 10110 : if (fockcommon%optstr.and.(fockcommon%ieigen/=0)) then
1323 : : #ifdef HAVE_OPENMP_OFFLOAD
1324 : : !$OMP TARGET EXIT DATA MAP(delete:strout) IF(gpu_option==ABI_GPU_OPENMP)
1325 : : #endif
1326 [ - + ]: 1870 : ABI_FREE(strout)
1327 : : end if
1328 : : #ifdef HAVE_OPENMP_OFFLOAD
1329 : : !$OMP TARGET EXIT DATA MAP(delete:grnhat_12) IF(gpu_option==ABI_GPU_OPENMP .and. ider==3)
1330 : : !$OMP TARGET EXIT DATA MAP(delete:rho12,gvnlxc) IF(gpu_option==ABI_GPU_OPENMP)
1331 : : !$OMP TARGET EXIT DATA MAP(delete:atom_ifftsph,atom_nfgd,atom_rfgd) IF(gpu_option==ABI_GPU_OPENMP)
1332 : : #endif
1333 [ - + ]: 10110 : ABI_FREE(grnhat_12)
1334 [ - + ]: 10110 : ABI_FREE(gvnlxc)
1335 [ - + ]: 10110 : ABI_FREE(grnhat12)
1336 [ - + ]: 10110 : ABI_FREE(rho12)
1337 [ - + ]: 10110 : ABI_FREE(atom_ifftsph)
1338 [ - + ]: 10110 : ABI_FREE(atom_nfgd)
1339 [ - + ]: 10110 : ABI_FREE(atom_rfgd)
1340 : : end if
1341 : :
1342 : 232015 : call timab(1528,2,tsec)
1343 : :
1344 : : end do ! jkpt
1345 : :
1346 : : ! ========================================================
1347 : : ! === After loop ===
1348 : : ! ========================================================
1349 : :
1350 : 29057 : call timab(1506,2,tsec) ; call timab(1507,1,tsec) ; call timab(1547,-1,tsec)
1351 : :
1352 [ + + ]: 29057 : if (fockcommon%usepaw==1) then
1353 [ + + + - ]: 2115 : if ((fockcommon%optfor).and.(fockcommon%ieigen/=0)) then
1354 : 1076 : call timab(1547,2,tsec) ; call timab(1548,-1,tsec)
1355 : 1076 : call xmpi_sum(forikpt,mpi_enreg%comm_hf,ier)
1356 : 1076 : call timab(1548,2,tsec) ; call timab(1547,-1,tsec)
1357 [ + + ]: 2152 : do idat=1,ndat
1358 [ + + ]: 4304 : do iatom=1,natom !Loop over atom
1359 : 2152 : ia=gs_ham%atindx(iatom)
1360 [ + - + + : 16140 : fockcommon%forces_ikpt(:,ia,fockcommon%ieigen+idat-1)=forikpt(:,iatom,idat)
+ + ]
1361 : : end do
1362 : : end do
1363 : : end if
1364 : : end if
1365 [ + + + - ]: 29057 : if(fockcommon%optstr.and.(fockcommon%ieigen/=0)) then
1366 : 2040 : call timab(1547,2,tsec) ; call timab(1548,-1,tsec)
1367 : 2040 : call xmpi_sum(fockcommon%stress_ikpt,mpi_enreg%comm_hf,ier)
1368 : 2040 : call timab(1548,2,tsec) ; call timab(1547,-1,tsec)
1369 : : end if
1370 : :
1371 [ + + ]: 29057 : if (.not.need_ghc) then
1372 : :
1373 : : ! ===============================
1374 : : ! === Deallocate local arrays ===
1375 : : ! ===============================
1376 : : #ifdef HAVE_OPENMP_OFFLOAD
1377 : : !$OMP TARGET EXIT DATA MAP(delete:cwavef_r,vqg,vlocpsi_r) IF(gpu_option==ABI_GPU_OPENMP)
1378 : : #endif
1379 [ - + ]: 2076 : ABI_FREE(cwavef_r)
1380 [ - + ]: 2076 : ABI_FREE(vlocpsi_r)
1381 [ - + ]: 2076 : ABI_FREE(dummytab)
1382 [ - + ]: 2076 : ABI_FREE(vqg)
1383 [ + + + - ]: 2076 : if(fockcommon%usepaw==1.or.fockcommon%optstr) then
1384 [ - + ]: 2076 : ABI_FREE(gboundf)
1385 : : end if
1386 [ + + ]: 2076 : if (fockcommon%usepaw==1) then
1387 [ + + + - ]: 238 : if ((fockcommon%optfor).and.(fockcommon%ieigen/=0)) then
1388 [ - + ]: 160 : ABI_FREE(forikpt)
1389 : : end if
1390 : : end if
1391 : : !*Restore gs_ham datastructure
1392 : :
1393 [ - + ]: 2076 : if (associated(gs_ham%ph3d_kp)) then
1394 [ # # # # : 0 : ABI_MALLOC(gs_ham%ph3d_kp,(2,gs_ham%npw_k,gs_ham%matblk))
# # # # #
# # # # #
# # ]
1395 : : end if
1396 : : call gs_ham%load_kprime(kpt_kp=gs_ham%kpt_k,istwf_kp=gs_ham%istwf_k,&
1397 : 2076 : & npw_kp=gs_ham%npw_k,kg_kp=gs_ham%kg_k,ffnl_kp=gs_ham%ffnl_k,ph3d_kp=gs_ham%ph3d_k)
1398 : :
1399 : : ! if (fockcommon%ieigen/=0) fockcommon%ieigen=0
1400 : :
1401 : : else
1402 : :
1403 : : ! *Restore gs_ham datastructure
1404 : :
1405 [ - + ]: 26981 : if (associated(gs_ham%ph3d_kp)) then
1406 [ # # # # : 0 : ABI_MALLOC(gs_ham%ph3d_kp,(2,gs_ham%npw_k,gs_ham%matblk))
# # # # #
# # # # #
# # ]
1407 : : end if
1408 : : call gs_ham%load_kprime(kpt_kp=gs_ham%kpt_k,istwf_kp=gs_ham%istwf_k,&
1409 : 26981 : & npw_kp=gs_ham%npw_k,kg_kp=gs_ham%kg_k,ffnl_kp=gs_ham%ffnl_k,ph3d_kp=gs_ham%ph3d_k)
1410 : :
1411 : : ! * Perform an FFT using fourwf to get ghc1 = FFT^-1(vlocpsi_r)
1412 : : #ifdef HAVE_OPENMP_OFFLOAD
1413 : : !$OMP TARGET UPDATE FROM(vlocpsi_r) IF(gpu_option==ABI_GPU_OPENMP)
1414 : : #endif
1415 [ - + + - : 134905 : ABI_MALLOC(psilocal,(cplex_fock*n4f,n5f,n6f*ndat))
- + + - -
+ + - + -
+ - - + -
+ ]
1416 [ + + ]: 53962 : do idat=1,ndat
1417 : : call fftpac(1,mpi_enreg,nspden_fock,cplex_fock*n1f,n2f,n3f,&
1418 : 53962 : cplex_fock*n4f,n5f,n6f,ngfft,vlocpsi_r(:,idat),psilocal(:,:,1+(idat-1)*n6f:idat*n6f),2)
1419 : : end do ! idat
1420 : :
1421 : 26981 : call timab(1515,2,tsec) ; call timab(1512,-1,tsec) ; call timab(1547,-2,tsec)
1422 : : call fourwf(0,rhodum0,rhodum,ghc1,psilocal,gboundf,gboundf,gs_ham%istwf_k,gs_ham%kg_k,gs_ham%kg_k,&
1423 : : & mgfftf,mpi_enreg,ndat,ngfftf,1,npw,n4f,n5f,n6f,3,tim_fourwf_fock_getghc,weight1,weight1,&
1424 [ + - - - : 26981 : & gpu_option=gs_ham%gpu_option)
- - - - +
- - - - -
- - - + -
+ - + -
+ ]
1425 : 26981 : call timab(1512,2,tsec) ; call timab(1515,-1,tsec) ; call timab(1547,-1,tsec)
1426 [ - + ]: 26981 : ABI_FREE(psilocal)
1427 : :
1428 [ + - ]: 26981 : if(gpu_option==ABI_GPU_DISABLED) then
1429 [ - + + + : 20086450 : ghc1=ghc1*sqrt(gs_ham%ucvol)+ghc2
+ + ]
1430 [ # # ]: 0 : else if(gpu_option==ABI_GPU_OPENMP) then
1431 : 0 : call abi_xaxpy(npw*ndat,cucvol,ghc1,1,ghc2,1,x_cplx=2,gpu_option=gpu_option)
1432 : : #ifdef HAVE_OPENMP_OFFLOAD
1433 : : call gpu_copy(ghc1,ghc2,int(2,c_size_t)*npw*ndat)
1434 : : !$OMP TARGET UPDATE FROM(ghc1)
1435 : : #endif
1436 : : end if
1437 : :
1438 : : ! * If the calculation is parallelized, perform an MPI_allreduce to sum all the contributions in the array ghc
1439 [ + + + + ]: 20059469 : ghc(:,:)=ghc(:,:)/mpi_enreg%nproc_hf + ghc1(:,:)
1440 : :
1441 : 26981 : call timab(1547,2,tsec) ; call timab(1548,-1,tsec)
1442 : 26981 : call xmpi_sum(ghc,mpi_enreg%comm_hf,ier)
1443 : 26981 : call timab(1548,2,tsec) ; call timab(1547,-1,tsec)
1444 : :
1445 : : ! ===============================
1446 : : ! === Deallocate local PAW arrays ===
1447 : : ! ===============================
1448 [ + + ]: 26981 : if (fockcommon%usepaw==1) then
1449 [ + + + - ]: 1877 : if ((fockcommon%optfor).and.(fockcommon%ieigen/=0)) then
1450 [ - + ]: 916 : ABI_FREE(forikpt)
1451 : : end if
1452 : : end if
1453 [ + + - + ]: 26981 : if(fockcommon%usepaw==1.or.fockcommon%optstr) then
1454 [ - + ]: 1877 : ABI_FREE(gboundf)
1455 : : end if
1456 : : ! ============================================
1457 : : ! === Calculate the contribution to energy ===
1458 : : ! ============================================
1459 : : ! * Only the contribution when cwavef=cgocc_bz are calculated, in order to cancel exactly the self-interaction
1460 : : ! * at each convergence step. (consistent definition with the definition of hartree energy)
1461 [ + - ]: 26981 : if (fockcommon%ieigen/=0) then
1462 [ + + ]: 53962 : do idat=1,ndat
1463 : : eigen=zero
1464 : : ! * Dot product of cwavef and ghc
1465 : : ! * inspired from the routine 54_spacepar/meanvalue_g but without the reference to parallelism and filtering
1466 [ + + ]: 26981 : if(gs_ham%istwf_k==2) then
1467 : 6 : eigen=half*cwavef(1,1+(idat-1)*npw)*ghc1(1,1+(idat-1)*npw)
1468 : : else
1469 : 26975 : eigen=cwavef(1,1+(idat-1)*npw)*ghc1(1,1+(idat-1)*npw)+cwavef(2,1+(idat-1)*npw)*ghc1(2,1+(idat-1)*npw)
1470 : : end if
1471 [ + + ]: 6677496 : do ipw=2,npw
1472 : 6677496 : eigen=eigen+cwavef(1,ipw+(idat-1)*npw)*ghc1(1,ipw+(idat-1)*npw)+cwavef(2,ipw+(idat-1)*npw)*ghc1(2,ipw+(idat-1)*npw)
1473 : : end do
1474 [ + + ]: 26981 : if(gs_ham%istwf_k>=2) eigen=two*eigen
1475 : 26981 : call timab(1547,2,tsec) ; call timab(1548,-1,tsec)
1476 : 26981 : call xmpi_sum(eigen,mpi_enreg%comm_hf,ier)
1477 : 26981 : call timab(1548,2,tsec) ; call timab(1547,-1,tsec)
1478 : 26981 : fockcommon%eigen_ikpt(fockcommon%ieigen+idat-1)= eigen
1479 [ - + ]: 80943 : if(fockcommon%use_ACE==0) fockcommon%ieigen = 0
1480 : : end do ! idat
1481 : : end if
1482 : :
1483 : : ! ===============================
1484 : : ! === Deallocate local arrays ===
1485 : : ! ===============================
1486 : : #ifdef HAVE_OPENMP_OFFLOAD
1487 : : !$OMP TARGET EXIT DATA MAP(delete:cwavef_r,ghc1,ghc2,vqg,vlocpsi_r) IF(gpu_option==ABI_GPU_OPENMP)
1488 : : #endif
1489 [ - + ]: 26981 : ABI_FREE(cwavef_r)
1490 [ - + ]: 26981 : ABI_FREE(ghc1)
1491 [ - + ]: 26981 : ABI_FREE(ghc2)
1492 [ - + ]: 26981 : ABI_FREE(vlocpsi_r)
1493 [ - + ]: 26981 : ABI_FREE(dummytab)
1494 [ - + ]: 26981 : ABI_FREE(vqg)
1495 : :
1496 : : endif
1497 : :
1498 : 29057 : call timab(1504,2,tsec) ; call timab(1507,-2,tsec) ; call timab(1515,-2,tsec) ; call timab(1547,-2,tsec)
1499 : :
1500 [ - + - + : 58114 : end subroutine fock_getghc
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + ]
1501 : : !!***
1502 : :
1503 : : !!****f* ABINIT/fock2ACE
1504 : : !! NAME
1505 : : !! fock2ACE
1506 : : !!
1507 : : !! FUNCTION
1508 : : !! Compute nonlocal contribution to the Fock part of the hamiltonian in the ACE formalism.
1509 : : !! optionally contribution to Fock forces
1510 : : !!
1511 : : !! INPUTS
1512 : : !! cg(2,mcg)=wavefunctions (may be read from disk file)
1513 : : !! cprj(natom,mcprj*usecprj)=<p_lmn|Cnk> coefficients for each WF |Cnk> and each NL proj |p_lmn>
1514 : : !! fock <type(fock_type)>= quantities to calculate Fock exact exchange
1515 : : !! istwfk(nkpt)=input option parameter that describes the storage of wfs
1516 : : !! kg(3,mpw*mkmem)=reduced coordinates (integers) of G vecs in basis
1517 : : !! kpt(3,nkpt)=k points in reduced coordinates
1518 : : !! mband=maximum number of bands
1519 : : !! mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol
1520 : : !! mcprj=size of projected wave-functions array (cprj) =nspinor*mband*mkmem*nsppol
1521 : : !! mgfft=maximum size of 1D FFTs
1522 : : !! mkmem=number of k points treated by this node.
1523 : : !! mpi_enreg=information about MPI parallelization
1524 : : !! mpsang=
1525 : : !! mpw= maximum number of plane waves
1526 : : !! my_natom=number of atoms treated by current processor
1527 : : !! natom=number of atoms in cell.
1528 : : !! nband(nkpt)=number of bands at each k point
1529 : : !! nfft=number of FFT grid points
1530 : : !! ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
1531 : : !! nkpt=number of k points in Brillouin zone
1532 : : !! nloalg(3)=governs the choice of the algorithm for non-local operator.
1533 : : !! npwarr(nkpt)=number of planewaves in basis and boundary at each k
1534 : : !! nspden=Number of spin Density components
1535 : : !! nspinor=number of spinorial components of the wavefunctions
1536 : : !! nsppol=1 for unpolarized, 2 for spin-polarized
1537 : : !! ntypat=number of types of atoms
1538 : : !! occ(mband*nkpt*nsppol)=occupation numbers for each band over all k points
1539 : : !! optfor=1 if computation of forces is required
1540 : : !! paw_ij(my_natom*usepaw) <type(paw_ij_type)>=paw arrays given on (i,j) channels
1541 : : !! pawtab(ntypat*usepaw) <type(pawtab_type)>=paw tabulated starting data
1542 : : !! ph1d(2,3*(2*mgfft+1)*natom)=one-dimensional structure factor information
1543 : : !! psps <type(pseudopotential_type)>=variables related to pseudopotentials
1544 : : !! rprimd(3,3)=dimensional primitive translations in real space (bohr)
1545 : : !! typat(natom)=type of each atom
1546 : : !! usecprj=1 if cprj datastructure has been allocated
1547 : : !! gpu_option= GPU implementation to use, i.e. cuda, openMP, ... (0=not using GPU)
1548 : : !! wtk(nkpt)=weight associated with each k point
1549 : : !! xred(3,natom)=reduced dimensionless atomic coordinates
1550 : : !! ylm(mpw*mkmem,mpsang*mpsang*useylm)= real spherical harmonics for each G and k point
1551 : : !!
1552 : : !! OUTPUT
1553 : : !!
1554 : : !! fock%fockACE(ikpt,isppol)%xi
1555 : : !! if optfor=1, fock%fock_common%forces
1556 : : !!
1557 : : !! SOURCE
1558 : :
1559 : 1240 : subroutine fock2ACE(cg,cprj,fock,istwfk,kg,kpt,mband,mcg,mcprj,mgfft,mkmem,mpi_enreg,mpsang,&
1560 : 1240 : & mpw,my_natom,natom,nband,nfft,ngfft,nkpt,nloalg,npwarr,nspden,nspinor,nsppol,&
1561 : 1240 : & ntypat,occ,optfor,paw_ij,pawtab,ph1d,psps,rprimd,typat,usecprj,gpu_option,wtk,xred,ylm)
1562 : :
1563 : : !Arguments ------------------------------------
1564 : : !scalars
1565 : : integer,intent(in) :: mband,mcg,mcprj,mgfft,mkmem,mpsang,mpw,my_natom,natom,nfft,nkpt
1566 : : integer,intent(in) :: nspden,nsppol,nspinor,ntypat,optfor
1567 : : integer,intent(in) :: usecprj,gpu_option
1568 : : type(MPI_type),intent(inout) :: mpi_enreg
1569 : : type(pseudopotential_type),intent(in) :: psps
1570 : : !arrays
1571 : : integer,intent(in) :: istwfk(nkpt),kg(3,mpw*mkmem),nband(nkpt*nsppol)
1572 : : integer,intent(in) :: ngfft(18),nloalg(3),npwarr(nkpt)
1573 : : integer,intent(in) :: typat(natom)
1574 : : real(dp),intent(in) :: cg(2,mcg)
1575 : : real(dp),intent(in) :: kpt(3,nkpt)
1576 : : real(dp),intent(in) :: occ(mband*nkpt*nsppol),ph1d(2,3*(2*mgfft+1)*natom)
1577 : : real(dp),intent(in) :: rprimd(3,3),wtk(nkpt),xred(3,natom)
1578 : : real(dp),intent(in) :: ylm(mpw*mkmem,mpsang*mpsang*psps%useylm)
1579 : : type(pawcprj_type),intent(inout) :: cprj(natom,mcprj*usecprj)
1580 : : type(paw_ij_type),intent(in) :: paw_ij(my_natom*psps%usepaw)
1581 : : type(pawtab_type),intent(in) :: pawtab(ntypat*psps%usepaw)
1582 : : type(fock_type),pointer, intent(inout) :: fock
1583 : : !Local variables-------------------------------
1584 : : !scalars
1585 : : integer :: bandpp,bdtot_index,dimffnl,iband,iband_cprj,iband_last,ibg,icg,ider
1586 : : integer :: ierr,info,idir,ikg,ikpt,ilm,ipw,isppol,istwf_k,kk,ll
1587 : : integer :: mband_cprj,me_distrb,my_ikpt,my_nspinor,nband_k,nband_cprj_k,nkpg
1588 : : integer :: npw_k,spaceComm
1589 : : integer :: use_ACE_old
1590 : : integer :: blocksize,iblock,jblock,iblocksize,jblocksize,nblockbd
1591 : : !integer, save :: counter=0
1592 : 1240 : type(gs_hamiltonian_type) :: gs_hamk
1593 : : logical :: compute_gbound
1594 : : type(fock_common_type),pointer :: fockcommon
1595 : : !arrays
1596 : 1240 : integer,allocatable :: kg_k(:,:)
1597 : : real(dp) :: kpoint(3),rmet(3,3),tsec(2)
1598 : 1240 : real(dp),allocatable :: bb(:,:,:),cwavef(:,:),cwavefk(:,:),ffnl_sav(:,:,:,:)
1599 : 1240 : real(dp),allocatable :: kpg_k(:,:),kpg_k_sav(:,:)
1600 : 1240 : real(dp),allocatable :: mkl(:,:,:),occblock(:),ph3d(:,:,:),ph3d_sav(:,:,:)
1601 : 1240 : real(dp),allocatable :: wi(:,:,:),weight(:),ylm_k(:,:),ylmgr_k(:,:,:)
1602 : 1240 : real(dp),allocatable,target :: ffnl(:,:,:,:)
1603 : : type(bandfft_kpt_type),pointer :: my_bandfft_kpt => null()
1604 : 1240 : type(pawcprj_type),target,allocatable :: cwaveprj(:,:)
1605 : : !*************************************************************************
1606 : :
1607 : 1240 : call timab(1560,1,tsec)
1608 : 1240 : call timab(1561,1,tsec)
1609 : :
1610 : : !DEBUG
1611 : : !if(counter>0)return
1612 : : !counter=counter+1
1613 : : !ENDDEBUG
1614 : :
1615 : : !Init mpicomm and me
1616 [ - + ]: 1240 : if(mpi_enreg%paral_kgb==1)then
1617 : 0 : spaceComm=mpi_enreg%comm_kpt
1618 : 0 : me_distrb=mpi_enreg%me_kpt
1619 : : else
1620 : : !* In case of HF calculation
1621 [ - + ]: 1240 : if (mpi_enreg%paral_hf==1) then
1622 : 0 : spaceComm=mpi_enreg%comm_kpt
1623 : 0 : me_distrb=mpi_enreg%me_kpt
1624 : : else
1625 : 1240 : spaceComm=mpi_enreg%comm_cell
1626 : 1240 : me_distrb=mpi_enreg%me_cell
1627 : : end if
1628 : : end if
1629 : :
1630 : : !Some initializations
1631 : 1240 : my_nspinor=max(1,nspinor/mpi_enreg%nproc_spinor)
1632 : 1240 : compute_gbound=.true.
1633 : 1240 : fockcommon => fock%fock_common
1634 : 1240 : use_ACE_old=fockcommon%use_ACE
1635 : 1240 : fockcommon%use_ACE=0
1636 : 1240 : fockcommon%e_fock0=zero
1637 : :
1638 : : !Initialize Hamiltonian (k- and spin-independent terms)
1639 : :
1640 : : call gs_hamk%init(psps,pawtab,nspinor,nsppol,nspden,natom,&
1641 : : & typat,xred,nfft,mgfft,ngfft,rprimd,nloalg,usecprj=usecprj,&
1642 : : & comm_atom=mpi_enreg%comm_atom,mpi_atmtab=mpi_enreg%my_atmtab,mpi_spintab=mpi_enreg%my_isppoltab,&
1643 : : & paw_ij=paw_ij,ph1d=ph1d,fock=fock,&
1644 : 1240 : & gpu_option=gpu_option)
1645 [ + + + + : 49600 : rmet = MATMUL(TRANSPOSE(rprimd),rprimd)
+ + ]
1646 : 1240 : fockcommon%use_ACE=use_ACE_old
1647 : :
1648 : : !need to reorder cprj=<p_lmn|Cnk> (from unsorted to atom-sorted)
1649 [ + + ]: 1240 : if (psps%usepaw==1) then
1650 : 206 : call pawcprj_reorder(cprj,gs_hamk%atindx)
1651 : : end if
1652 : :
1653 : : !LOOP OVER SPINS
1654 : 1240 : bdtot_index=0;ibg=0;icg=0
1655 : 1240 : call timab(1561,2,tsec) ; call timab(1562,-1,tsec)
1656 : :
1657 [ + + ]: 2615 : do isppol=1,nsppol
1658 : 1375 : fockcommon%isppol=isppol
1659 : : ! Continue to initialize the Hamiltonian (PAW DIJ coefficients)
1660 : 1375 : call gs_hamk%load_spin(isppol,with_nonlocal=.true.)
1661 : :
1662 : : ! Loop over k points
1663 : 1375 : ikg=0
1664 [ + + ]: 9919 : do ikpt=1,nkpt
1665 : 7304 : fockcommon%ikpt=ikpt
1666 : 7304 : nband_k=nband(ikpt+(isppol-1)*nkpt)
1667 : 7304 : npw_k=npwarr(ikpt)
1668 [ + + ]: 29216 : kpoint(:)=kpt(:,ikpt)
1669 : 7304 : istwf_k=istwfk(ikpt)
1670 [ + + ]: 7304 : if(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,isppol,me_distrb)) then
1671 : 1720 : bdtot_index=bdtot_index+nband_k
1672 : 1720 : cycle
1673 : : end if
1674 : :
1675 : : ! Parallelism over FFT and/or bands: define sizes and tabs
1676 [ - + ]: 5584 : if (mpi_enreg%paral_kgb==1) then
1677 : 0 : my_ikpt=mpi_enreg%my_kpttab(ikpt)
1678 : 0 : nblockbd=nband_k/(mpi_enreg%nproc_band*mpi_enreg%bandpp)
1679 : 0 : bandpp=mpi_enreg%bandpp
1680 : 0 : my_bandfft_kpt => bandfft_kpt(my_ikpt)
1681 : : else
1682 : 5584 : my_ikpt=ikpt
1683 : 5584 : bandpp=mpi_enreg%bandpp
1684 : 5584 : nblockbd=nband_k/bandpp
1685 : : end if
1686 : 5584 : blocksize=nband_k/nblockbd
1687 : 5584 : mband_cprj=mband/mpi_enreg%nproc_band
1688 : 5584 : nband_cprj_k=nband_k/mpi_enreg%nproc_band
1689 : :
1690 [ - + + - : 16752 : ABI_MALLOC(cwavef,(2,npw_k*my_nspinor*blocksize))
+ - + - -
+ - + -
+ ]
1691 [ + + ]: 5584 : if (psps%usepaw==1) then
1692 [ - + + - : 2954 : ABI_MALLOC(cwaveprj,(natom,my_nspinor*bandpp))
- + + - +
- + - - +
- + - + +
+ + + ]
1693 : 422 : call pawcprj_alloc(cwaveprj,0,gs_hamk%dimcprj)
1694 : : else
1695 [ - + - + ]: 5162 : ABI_MALLOC(cwaveprj,(0,0))
1696 : : end if
1697 : :
1698 [ - + + - : 16752 : ABI_MALLOC(kg_k,(3,mpw))
+ - + - -
+ - + -
+ ]
1699 : : !$OMP PARALLEL DO
1700 [ + + ]: 1241619 : do ipw=1,npw_k
1701 [ + + ]: 4949724 : kg_k(:,ipw)=kg(:,ipw+ikg)
1702 : : end do
1703 : :
1704 [ - + + - : 22336 : ABI_MALLOC(ylm_k,(npw_k,mpsang*mpsang*psps%useylm))
+ + + - +
- + + - +
- + - + ]
1705 [ - + - + ]: 5584 : ABI_MALLOC(ylmgr_k,(0,0,0))
1706 [ + + ]: 5584 : if (psps%useylm==1) then
1707 : : !$OMP PARALLEL DO COLLAPSE(2)
1708 [ + + ]: 2110 : do ilm=1,mpsang*mpsang
1709 [ + + ]: 1282866 : do ipw=1,npw_k
1710 : 1282444 : ylm_k(ipw,ilm)=ylm(ipw+ikg,ilm)
1711 : : end do
1712 : : end do
1713 : : end if
1714 : :
1715 : : ! Compute (k+G) vectors
1716 : 5584 : nkpg=3*nloalg(3)
1717 [ - + + - : 22336 : ABI_MALLOC(kpg_k,(npw_k,nkpg))
+ + + - +
- + + - +
- + - + ]
1718 [ + + ]: 5584 : if (nkpg>0) then
1719 : 422 : call mkkpg(kg_k,kpg_k,kpoint,nkpg,npw_k)
1720 : : end if
1721 : :
1722 : :
1723 : : ! Compute nonlocal form factors ffnl at all (k+G)
1724 : 5584 : ider=0;idir=0;dimffnl=1
1725 [ - + + - : 27920 : ABI_MALLOC(ffnl,(npw_k,dimffnl,psps%lmnmax,ntypat))
- + + - -
+ + - + -
+ - - + -
+ - + ]
1726 : : call mkffnl(psps%dimekb,dimffnl,psps%ekb,ffnl,psps%ffspl,gs_hamk%gmet,gs_hamk%gprimd,&
1727 : : & ider,idir,psps%indlmn,kg_k,kpg_k,kpoint,psps%lmnmax,psps%lnmax,psps%mpsang,psps%mqgrid_ff,&
1728 : 5584 : & nkpg,npw_k,ntypat,psps%pspso,psps%qgrid_ff,rmet,psps%usepaw,psps%useylm,ylm_k,ylmgr_k)
1729 : :
1730 : : ! Load k-dependent part in the Hamiltonian datastructure
1731 : : ! - Compute 3D phase factors
1732 : : ! - Prepare various tabs in case of band-FFT parallelism
1733 : : ! - Load k-dependent quantities in the Hamiltonian
1734 : :
1735 [ - + + - : 22336 : ABI_MALLOC(ph3d,(2,npw_k,gs_hamk%matblk))
- + + - +
- + - - +
- + - + ]
1736 : : call gs_hamk%load_k(kpt_k=kpoint,istwf_k=istwf_k,npw_k=npw_k,&
1737 : 5584 : & kg_k=kg_k,kpg_k=kpg_k,ffnl_k=ffnl,ph3d_k=ph3d,compute_gbound=compute_gbound,compute_ph3d=.true.)
1738 : :
1739 : : ! Load band-FFT tabs (transposed k-dependent arrays)
1740 [ - + ]: 5584 : if (mpi_enreg%paral_kgb==1) then
1741 : 0 : call bandfft_kpt_savetabs(my_bandfft_kpt,ffnl=ffnl_sav,ph3d=ph3d_sav,kpg=kpg_k_sav)
1742 : 0 : call prep_bandfft_tabs(gs_hamk,ikpt,mkmem,mpi_enreg)
1743 : : call gs_hamk%load_k(npw_fft_k=my_bandfft_kpt%ndatarecv, &
1744 : : & kg_k =my_bandfft_kpt%kg_k_gather, &
1745 : : & kpg_k =my_bandfft_kpt%kpg_k_gather, &
1746 : : ffnl_k =my_bandfft_kpt%ffnl_gather, &
1747 : 0 : ph3d_k =my_bandfft_kpt%ph3d_gather,compute_gbound=compute_gbound)
1748 : : end if
1749 : :
1750 : : ! The following is now wrong. In sequential, nblockbd=nband_k/bandpp
1751 : : ! blocksize= bandpp (JB 2016/04/16)
1752 : : ! Note that in sequential mode iblock=iband, nblockbd=nband_k and blocksize=1
1753 : : !
1754 [ - + + - : 16752 : ABI_MALLOC(occblock,(blocksize))
+ - + - -
+ - + ]
1755 [ - + + - : 11168 : ABI_MALLOC(weight,(blocksize))
- + - + ]
1756 [ + + + + ]: 16752 : occblock=zero;weight=zero
1757 : :
1758 [ + + ]: 5584 : if (fockcommon%optfor) then
1759 [ + + + + : 13345 : fockcommon%forces_ikpt=zero
+ + ]
1760 : : end if
1761 : :
1762 [ - + + - : 22336 : ABI_MALLOC(wi,(2,npw_k*my_nspinor*blocksize,nblockbd))
- + + - +
- + - - +
- + ]
1763 [ + + + + : 20065053 : wi=zero
+ + ]
1764 [ - + + - : 22336 : ABI_MALLOC(mkl,(2,nband_k,nband_k))
- + + - +
- + - - +
- + - + ]
1765 [ + + + + : 541824 : mkl=zero
+ + ]
1766 : : ! Calculate all the Wi for the current k-point
1767 : :
1768 [ + + ]: 32565 : do iblock=1,nblockbd
1769 : :
1770 : 26981 : iband=(iblock-1)*blocksize+1;iband_last=min(iband+blocksize-1,nband_k)
1771 : 26981 : iband_cprj=(iblock-1)*bandpp+1
1772 [ + - ]: 26981 : if(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,iband,iband_last,isppol,me_distrb)) cycle
1773 : :
1774 : : ! Select occupied bandsddk
1775 [ + + ]: 53962 : occblock(:)=occ(1+(iblock-1)*blocksize+bdtot_index:iblock*blocksize+bdtot_index)
1776 [ + + ]: 53962 : weight(:)=wtk(ikpt)*occblock(:)
1777 : :
1778 : : ! Load contribution from n,k
1779 : : cwavef(:,1:npw_k*my_nspinor*blocksize)=&
1780 [ + + + + ]: 20059469 : & cg(:,1+(iblock-1)*npw_k*my_nspinor*blocksize+icg:iblock*npw_k*my_nspinor*blocksize+icg)
1781 [ + + ]: 26981 : if (psps%usepaw==1) then
1782 : : call pawcprj_get(gs_hamk%atindx1,cwaveprj,cprj,natom,iband_cprj,ibg,ikpt,0,isppol,&
1783 : : & mband_cprj,mkmem,natom,bandpp,nband_cprj_k,my_nspinor,nsppol,0,&
1784 : 1877 : & mpicomm=mpi_enreg%comm_kpt,proc_distrb=mpi_enreg%proc_distrb)
1785 : : end if
1786 : :
1787 [ - + ]: 26981 : if (mpi_enreg%paral_kgb==1) then
1788 : 0 : ABI_BUG('Paral_kgb is not yet implemented for fock calculations')
1789 : : end if
1790 : :
1791 : 26981 : fockcommon%ieigen=(iblock-1)*blocksize+1
1792 : 26981 : fockcommon%iband=(iblock-1)*blocksize+1
1793 : 26981 : call timab(1562,2,tsec) ; call timab(1563,-1,tsec)
1794 : : call fock_getghc(cwavef,&
1795 : : cwaveprj,&
1796 : : wi(:,:,iblock),&
1797 : 26981 : gs_hamk,mpi_enreg,blocksize)
1798 : :
1799 : :
1800 : 26981 : call timab(1563,2,tsec) ; call timab(1562,-1,tsec)
1801 [ + + ]: 59546 : do iblocksize=1,blocksize
1802 : 26981 : mkl(1,fockcommon%ieigen+iblocksize-1,fockcommon%ieigen+iblocksize-1)=fockcommon%eigen_ikpt(fockcommon%ieigen+iblocksize-1)
1803 : 26981 : fockcommon%e_fock0=fockcommon%e_fock0+half*weight(iblocksize)*fockcommon%eigen_ikpt(fockcommon%ieigen+iblocksize-1)
1804 [ + + ]: 53962 : if (fockcommon%optfor) then
1805 [ + + + + ]: 12996 : fockcommon%forces(:,:)=fockcommon%forces(:,:)+weight(iblocksize)*fockcommon%forces_ikpt(:,:,fockcommon%ieigen+iblocksize-1)
1806 : : end if
1807 : : end do
1808 : :
1809 : :
1810 : : end do ! End of loop on block of bands
1811 : :
1812 : : ! Calculate Mkl for the current k-point
1813 [ - + + - : 16752 : ABI_MALLOC(cwavefk,(2,npw_k*my_nspinor))
+ - + - -
+ - + ]
1814 [ + + ]: 32565 : do iblock=1,nblockbd
1815 : : cwavef(:,1:npw_k*my_nspinor*blocksize)=&
1816 [ + + + + ]: 20059469 : & cg(:,1+(iblock-1)*npw_k*my_nspinor*blocksize+icg:iblock*npw_k*my_nspinor*blocksize+icg)
1817 [ + + ]: 59546 : do iblocksize=1,blocksize
1818 : 26981 : kk=(iblock-1)*blocksize+iblocksize
1819 [ + + + + ]: 20059469 : cwavefk(:,:)=cwavef(:,1+(iblocksize-1)*npw_k*my_nspinor:iblocksize*npw_k*my_nspinor)
1820 [ + + ]: 152329 : do jblock=1,iblock
1821 [ + + ]: 223715 : do jblocksize=1,blocksize
1822 : 98367 : ll=(jblock-1)*blocksize+jblocksize
1823 [ + + ]: 196734 : if (ll<kk) then
1824 : : call dotprod_g(mkl(1,kk,ll),mkl(2,kk,ll),gs_hamk%istwf_k,npw_k,2,wi(:,1+(jblocksize-1)*npw_k*my_nspinor:&
1825 : 71386 : & jblocksize*npw_k*my_nspinor,jblock),cwavefk,mpi_enreg%me_g0,mpi_enreg%comm_fft)
1826 : : end if
1827 : : end do
1828 : : end do
1829 : : end do
1830 : : end do ! End of loop on block of bands
1831 : :
1832 [ - + ]: 5584 : ABI_FREE(cwavefk)
1833 [ + + + + : 541824 : mkl=-mkl
+ + ]
1834 : :
1835 : : ! Cholesky factorisation of -mkl=Lx(trans(L)*. On output mkl=L
1836 : 5584 : call abi_xpotrf("L",nband_k,mkl,nband_k,info,x_cplx=2)
1837 : :
1838 : : ! calculate trans(L-1)
1839 [ - + + - : 16752 : ABI_MALLOC(bb,(2,nband_k,nband_k))
- + + - -
+ - + ]
1840 [ + + + + : 541824 : bb=zero
+ + ]
1841 [ + + ]: 32565 : do kk=1,nband_k
1842 : 32565 : bb(1,kk,kk)=one
1843 : : end do
1844 : 5584 : call ztrtrs("L","T","N",nband_k,nband_k,mkl,nband_k,bb,nband_k,info)
1845 [ + + + + : 20065053 : fock%fockACE(ikpt,isppol)%xi=zero
+ + ]
1846 : :
1847 : : ! Calculate ksi
1848 [ + + ]: 32565 : do kk=1,nband_k
1849 [ + + ]: 202318 : do jblock=1,nblockbd
1850 [ + + ]: 366487 : do jblocksize=1,blocksize
1851 : 169753 : ll=(jblock-1)*blocksize+jblocksize
1852 : : fock%fockACE(ikpt,isppol)%xi(1,:,kk)=fock%fockACE(ikpt,isppol)%xi(1,:,kk)+bb(1,ll,kk)*wi(1,1+(jblocksize-1)*&
1853 : : & npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock)-&
1854 [ + + ]: 42673461 : & bb(2,ll,kk)*wi(2,1+(jblocksize-1)*npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock)
1855 : : fock%fockACE(ikpt,isppol)%xi(2,:,kk)=fock%fockACE(ikpt,isppol)%xi(2,:,kk)+bb(1,ll,kk)*wi(2,1+(jblocksize-1)*&
1856 : : npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock)+&
1857 [ + + ]: 42843214 : & bb(2,ll,kk)*wi(1,1+(jblocksize-1)*npw_k*my_nspinor:jblocksize*npw_k*my_nspinor,jblock)
1858 : : end do
1859 : : end do
1860 : : end do
1861 : :
1862 : : ! DEBUG
1863 : : ! fock%fockACE(ikpt,isppol)%xi=zero
1864 : : ! ENDDEBUG
1865 : :
1866 [ - + ]: 5584 : ABI_FREE(wi)
1867 [ - + ]: 5584 : ABI_FREE(mkl)
1868 : :
1869 : : ! Restore the bandfft tabs
1870 [ - + ]: 5584 : if (mpi_enreg%paral_kgb==1) then
1871 : 0 : call bandfft_kpt_restoretabs(my_bandfft_kpt,ffnl=ffnl_sav,ph3d=ph3d_sav,kpg=kpg_k_sav)
1872 : : end if
1873 : :
1874 : : ! Increment indices
1875 : 5584 : bdtot_index=bdtot_index+nband_k
1876 [ + - ]: 5584 : if (mkmem/=0) then
1877 : 5584 : ibg=ibg+my_nspinor*nband_cprj_k
1878 : 5584 : icg=icg+npw_k*my_nspinor*nband_k
1879 : 5584 : ikg=ikg+npw_k
1880 : : end if
1881 : :
1882 [ + + ]: 5584 : if (psps%usepaw==1) then
1883 : 422 : call pawcprj_free(cwaveprj)
1884 : : end if
1885 [ + - + + : 6428 : ABI_FREE(cwaveprj)
- + - + -
+ ]
1886 [ - + ]: 5584 : ABI_FREE(cwavef)
1887 [ - + ]: 5584 : ABI_FREE(bb)
1888 [ - + ]: 5584 : ABI_FREE(occblock)
1889 [ - + ]: 5584 : ABI_FREE(weight)
1890 [ - + ]: 5584 : ABI_FREE(ffnl)
1891 [ - + ]: 5584 : ABI_FREE(kg_k)
1892 [ - + ]: 5584 : ABI_FREE(kpg_k)
1893 [ - + ]: 5584 : ABI_FREE(ylm_k)
1894 [ - + ]: 5584 : ABI_FREE(ylmgr_k)
1895 [ - + ]: 6959 : ABI_FREE(ph3d)
1896 : : end do ! End k point loop
1897 : : end do ! End loop over spins
1898 : :
1899 : 1240 : call timab(1562,2,tsec)
1900 : 1240 : call timab(1565,1,tsec)
1901 : :
1902 : : !Parallel case: accumulate (n,k) contributions
1903 : : if (xmpi_paral==1) then
1904 : 1240 : call xmpi_sum(fockcommon%e_fock0,spaceComm,ierr)
1905 : : ! Forces
1906 [ + + ]: 1240 : if (optfor==1) then
1907 [ + + ]: 220 : if (psps%usepaw==1) then
1908 : 144 : call xmpi_sum(fockcommon%forces,spaceComm,ierr)
1909 : : end if
1910 : : end if
1911 : : end if
1912 : :
1913 : : !need to reorder cprj=<p_lmn|Cnk> (from atom-sorted to unsorted)
1914 [ + + ]: 1240 : if (psps%usepaw==1) then
1915 : 206 : call pawcprj_reorder(cprj,gs_hamk%atindx1)
1916 : : end if
1917 : : !Deallocate temporary space
1918 : 1240 : call gs_hamk%free()
1919 : :
1920 : 1240 : call timab(1565,2,tsec)
1921 : 1240 : call timab(1560,2,tsec)
1922 : :
1923 [ - + - + : 2480 : end subroutine fock2ACE
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - - - -
- - - + ]
1924 : : !!***
1925 : :
1926 : : !!****f* ABINIT/fock_ACE_getghc
1927 : : !! NAME
1928 : : !! fock_ACE_getghc
1929 : : !!
1930 : : !! FUNCTION
1931 : : !! Compute the matrix elements <G|Vx|psi> of the Fock operator in the ACE context.
1932 : : !!
1933 : : !! INPUTS
1934 : : !! cwavef(2,npw*nspinor*ndat)= planewave coefficients of wavefunctions on which Fock operator is applied.
1935 : : !! gs_ham <type(gs_hamiltonian_type)>=all data for the Hamiltonian to be applied
1936 : : !! mpi_enreg= information about MPI parallelization
1937 : : !!
1938 : : !! SIDE EFFECTS
1939 : : !! ghc(2,npw*ndat)= matrix elements <G|H|C> or <G|H-lambda.S|C> (if sij_opt>=0 or =-1 in getghc)
1940 : : !! contains the fock exchange term for cwavef at the end.
1941 : : !!
1942 : : !! NOTES
1943 : : !! The current version assumes that :
1944 : : !! * nspinor = 1
1945 : : !! * no "my_nspinor"
1946 : : !! * no restriction to the value of istwfk_bz (but must be tested in all case)
1947 : : !! * all the data for the occupied states (cgocc_bz) are the same as those for the current states (cg)
1948 : : !!
1949 : : !! SOURCE
1950 : :
1951 [ - + - + ]: 142507 : subroutine fock_ACE_getghc(cwavef,ghc,gs_ham,mpi_enreg,ndat,gpu_option)
1952 : :
1953 : : !Arguments ------------------------------------
1954 : : ! Scalars
1955 : : integer :: ndat
1956 : : integer,optional :: gpu_option
1957 : : type(MPI_type),intent(in) :: mpi_enreg
1958 : : type(gs_hamiltonian_type),target,intent(inout) :: gs_ham
1959 : : ! Arrays
1960 : : real(dp),target,intent(inout) :: cwavef(:,:)!,ghc(2,gs_ham%npw_k*ndat)
1961 : : real(dp),target,intent(inout) :: ghc(:,:)
1962 : :
1963 : : !Local variables-------------------------------
1964 : : ! Scalars
1965 : : complex(dp), parameter :: cminusone = (-1._dp,0._dp)
1966 : : integer :: iband,ikpt,ipw,my_nspinor,nband_k,npw,idat,gpu_option_
1967 : : real(dp) :: eigen
1968 : : type(fock_common_type),pointer :: fockcommon
1969 : : ! Arrays
1970 : : real(dp) :: tsec(2)
1971 : 142507 : real(dp), target, allocatable :: mat(:,:,:),ghc1(:,:)
1972 : 142507 : real(dp), allocatable :: vdotr(:),vdoti(:)
1973 : 142507 : real(dp), contiguous, pointer :: xi(:,:,:)
1974 : : ! *************************************************************************
1975 : :
1976 : 142507 : call timab(1580,1,tsec)
1977 : :
1978 [ - + ]: 142507 : ABI_CHECK(associated(gs_ham%fockcommon),"fock must be associated!")
1979 : 142507 : fockcommon => gs_ham%fockcommon
1980 : :
1981 [ - + ]: 142507 : ABI_CHECK(gs_ham%nspinor==1,"only allowed for nspinor=1!")
1982 [ - + ]: 142507 : ABI_CHECK(gs_ham%npw_k==gs_ham%npw_kp,"only allowed for npw_k=npw_kp (ground state)!")
1983 : :
1984 : 142507 : ikpt=fockcommon%ikpt
1985 : 142507 : npw=gs_ham%npw_k
1986 : 142507 : nband_k=fockcommon%nband(ikpt)
1987 : 142507 : my_nspinor=max(1,gs_ham%nspinor/mpi_enreg%nproc_spinor)
1988 [ - + ]: 142507 : gpu_option_=ABI_GPU_DISABLED; if(present(gpu_option)) gpu_option_ = gpu_option
1989 : : !*Initialization of the array ghc1
1990 : : !*ghc1 will contain the exact exchange contribution to the Hamiltonian
1991 [ - + + - : 427521 : ABI_MALLOC(ghc1,(2,npw*my_nspinor*ndat))
+ - + - -
+ - + ]
1992 [ + + + + ]: 107530849 : ghc1=zero
1993 : :
1994 : 142507 : xi => gs_ham%fockACE_k%xi(:,:,:)
1995 : :
1996 [ + - ]: 142507 : if(gpu_option_==ABI_GPU_DISABLED) then
1997 [ + + ]: 142507 : if(gs_ham%istwf_k==1) then
1998 [ - + + - : 569484 : ABI_MALLOC(mat,(2,nband_k,ndat))
- + + - +
- + - - +
- + ]
1999 : : call abi_zgemm_2r('C', 'N', nband_k, ndat, npw, cone, &
2000 : : xi, npw, &
2001 : : cwavef, npw, &
2002 : : czero, &
2003 [ - + - - : 142371 : mat, nband_k)
- - - - -
+ ]
2004 : : call abi_zgemm_2r('N', 'N', npw, ndat, nband_k, cminusone, &
2005 : : xi, npw, &
2006 : : mat, nband_k, &
2007 : : czero, &
2008 : 142371 : ghc1, npw)
2009 [ - + ]: 142371 : ABI_FREE(mat)
2010 : : else
2011 [ - + + - : 408 : ABI_MALLOC(vdotr,(nband_k))
+ - + - -
+ - + ]
2012 [ - + + - : 272 : ABI_MALLOC(vdoti,(nband_k))
- + + - ]
2013 [ + + ]: 272 : do idat=1,ndat
2014 : : call dotprod_g_batch_half(vdotr,vdoti,gs_ham%istwf_k,npw*my_nspinor,nband_k,2,&
2015 [ - + - - : 136 : & cwavef(:,1+(idat-1)*npw:idat*npw),xi(:,:,:),mpi_enreg%me_g0,mpi_enreg%comm_fft)
- - - - -
+ ]
2016 : :
2017 [ + + ]: 408 : do iband=1, nband_k
2018 [ + - + + : 820624 : ghc1(1,1+(idat-1)*npw:idat*npw)=ghc1(1,1+(idat-1)*npw:idat*npw)-vdotr(iband)*xi(1,:,iband)
+ + ]
2019 [ + + + + ]: 820760 : ghc1(2,1+(idat-1)*npw:idat*npw)=ghc1(2,1+(idat-1)*npw:idat*npw)-vdotr(iband)*xi(2,:,iband)
2020 : : end do
2021 : : end do
2022 [ - + ]: 136 : ABI_FREE(vdotr)
2023 [ - + ]: 136 : ABI_FREE(vdoti)
2024 : : end if
2025 : :
2026 : : !* If the calculation is parallelized, perform an MPI_allreduce to sum all the contributions in the array ghc
2027 : : ! ghc(:,:)=ghc(:,:)/mpi_enreg%nproc_spkpt + ghc1(:,:)
2028 [ + - + + : 215061698 : ghc(:,:)=ghc(:,:) + ghc1(:,:)
+ + + + +
+ ]
2029 : :
2030 : : ! call xmpi_sum(ghc,mpi_enreg%comm_kpt,ier)
2031 : :
2032 : : else if(gpu_option_==ABI_GPU_OPENMP) then
2033 : : #ifdef HAVE_OPENMP_OFFLOAD
2034 : : !$OMP TARGET ENTER DATA MAP(alloc:xi,ghc1)
2035 : : if(gs_ham%istwf_k==1) then
2036 : : ABI_MALLOC(mat,(2,nband_k,ndat))
2037 : : !$OMP TARGET ENTER DATA MAP(alloc:mat)
2038 : : !$OMP TARGET UPDATE TO(xi)
2039 : : call abi_zgemm_2r('C', 'N', nband_k, ndat, npw, cone, &
2040 : : & xi, npw, &
2041 : : & cwavef, npw, &
2042 : : & czero, &
2043 : : & mat, nband_k, &
2044 : : & gpu_option=gpu_option_)
2045 : : call abi_zgemm_2r('N', 'N', npw, ndat, nband_k, cminusone, &
2046 : : & xi, npw, &
2047 : : & mat, nband_k, &
2048 : : & czero, &
2049 : : & ghc1, npw, &
2050 : : & gpu_option=gpu_option_)
2051 : :
2052 : : !$OMP TARGET EXIT DATA MAP(delete:mat)
2053 : : ABI_FREE(mat)
2054 : : else
2055 : : ABI_MALLOC(vdotr,(nband_k))
2056 : : ABI_MALLOC(vdoti,(nband_k))
2057 : : !$OMP TARGET ENTER DATA MAP(alloc:vdotr,vdoti)
2058 : : do idat=1,ndat
2059 : : call dotprod_g_batch_half(vdotr,vdoti,gs_ham%istwf_k,npw*my_nspinor,nband_k,2,&
2060 : : & cwavef(:,1+(idat-1)*npw:idat*npw),xi(:,:,:),mpi_enreg%me_g0,mpi_enreg%comm_fft,&
2061 : : & gpu_option=gpu_option_)
2062 : :
2063 : : !$OMP TARGET TEAMS DISTRIBUTE PRIVATE(iband) MAP(ghc1,xi,vdotr)
2064 : : do iband=1, nband_k
2065 : : !$OMP PARALLEL DO PRIVATE(ipw)
2066 : : do ipw=1,npw
2067 : : ghc1(1,ipw+(idat-1)*npw)=ghc1(1,ipw+(idat-1)*npw)-vdotr(iband)*xi(1,ipw,iband)
2068 : : ghc1(2,ipw+(idat-1)*npw)=ghc1(2,ipw+(idat-1)*npw)-vdotr(iband)*xi(2,ipw,iband)
2069 : : end do
2070 : : end do
2071 : : end do
2072 : : !$OMP TARGET EXIT DATA MAP(delete:vdotr,vdoti)
2073 : : ABI_FREE(vdotr)
2074 : : ABI_FREE(vdoti)
2075 : : end if
2076 : :
2077 : : !* If the calculation is parallelized, perform an MPI_allreduce to sum all the contributions in the array ghc
2078 : : ! ghc(:,:)=ghc(:,:)/mpi_enreg%nproc_spkpt + ghc1(:,:)
2079 : :
2080 : : call abi_xaxpy(npw*ndat,cone,ghc1,1,ghc,1,x_cplx=2,gpu_option=gpu_option)
2081 : :
2082 : : ! call xmpi_sum(ghc,mpi_enreg%comm_kpt,ier)
2083 : : !$OMP TARGET UPDATE FROM(ghc1)
2084 : :
2085 : : !$OMP TARGET EXIT DATA MAP(delete:xi,ghc1)
2086 : : #endif
2087 : : end if
2088 : :
2089 : :
2090 : : ! ============================================
2091 : : ! === Calculate the contribution to energy ===
2092 : : ! ============================================
2093 : : !* Only the contribution when cwavef=cgocc_bz are calculated, in order to cancel exactly the self-interaction
2094 : : !* at each convergence step. (consistent definition with the definition of hartree energy)
2095 [ + + ]: 142507 : if (fockcommon%ieigen/=0) then
2096 [ + + ]: 68638 : do idat=1,ndat
2097 : 34319 : eigen=zero
2098 : : ! * Dot product of cwavef and ghc
2099 : : ! * inspired from the routine 54_spacepar/meanvalue_g but without the reference to parallelism and filtering
2100 [ + + ]: 34319 : if(gs_ham%istwf_k==2) then
2101 : 29 : eigen=half*cwavef(1,1+(idat-1)*npw)*ghc1(1,1+(idat-1)*npw)
2102 : : else
2103 : 34290 : eigen=cwavef(1,1+(idat-1)*npw)*ghc1(1,1+(idat-1)*npw)+cwavef(2,1+(idat-1)*npw)*ghc1(2,1+(idat-1)*npw)
2104 : : end if
2105 [ + + ]: 8422706 : do ipw=2,npw
2106 : 8422706 : eigen=eigen+cwavef(1,ipw+(idat-1)*npw)*ghc1(1,ipw+(idat-1)*npw)+cwavef(2,ipw+(idat-1)*npw)*ghc1(2,ipw+(idat-1)*npw)
2107 : : end do
2108 [ + + ]: 34319 : if(gs_ham%istwf_k>=2) eigen=two*eigen
2109 : : ! call xmpi_sum(eigen,mpi_enreg%comm_kpt,ier)
2110 : 34319 : fockcommon%eigen_ikpt(fockcommon%ieigen+idat-1)= eigen
2111 : 68638 : fockcommon%ieigen = 0
2112 : : end do ! idat
2113 : : end if
2114 : :
2115 : : ! ===============================
2116 : : ! === Deallocate local arrays ===
2117 : : ! ===============================
2118 : :
2119 [ - + ]: 142507 : ABI_FREE(ghc1)
2120 : 142507 : call timab(1580,2,tsec)
2121 : :
2122 : 142507 : end subroutine fock_ACE_getghc
2123 : : !!***
2124 : :
2125 : : end module m_fock_getghc
2126 : : !!***
|
