Branch data Line data Source code
1 : : !!****f* ABINIT/m_chebfiwf
2 : : !! NAME
3 : : !! m_chebfiwf
4 : : !!
5 : : !! FUNCTION
6 : : !! This module contains a routine updating the whole wave functions at a given k-point,
7 : : !! using the Chebyshev filtering method (2021 implementation using xG abstraction layer)
8 : : !! for a given spin-polarization, from a fixed hamiltonian
9 : : !! but might also simply compute eigenvectors and eigenvalues at this k point.
10 : : !! it will also update the matrix elements of the hamiltonian.
11 : : !!
12 : : !! COPYRIGHT
13 : : !! Copyright (C) 2018-2026 ABINIT group (BS, IML)
14 : : !! This file is distributed under the terms of the
15 : : !! gnu general public license, see ~abinit/COPYING
16 : : !! or http://www.gnu.org/copyleft/gpl.txt .
17 : : !! for the initials of contributors, see ~abinit/doc/developers/contributors.txt .
18 : : !!
19 : : !! SOURCE
20 : :
21 : : #if defined HAVE_CONFIG_H
22 : : #include "config.h"
23 : : #endif
24 : :
25 : : #include "abi_common.h"
26 : :
27 : : ! nvtx related macro definition
28 : : #include "nvtx_macros.h"
29 : :
30 : : module m_chebfiwf
31 : :
32 : : use defs_abitypes
33 : : use defs_basis
34 : : use m_abicore
35 : : use m_errors
36 : : use m_fstrings
37 : : use m_time
38 : :
39 : : use m_chebfi
40 : : use m_chebfi2
41 : : use m_invovl
42 : :
43 : : use m_cgtools, only : dotprod_g
44 : : use m_dtset, only : dataset_type
45 : :
46 : : use m_hamiltonian, only : gs_hamiltonian_type
47 : : use m_pawcprj, only : pawcprj_type
48 : : use m_nonlop, only : nonlop
49 : : use m_prep_kgb, only : prep_getghc, prep_nonlop
50 : : use m_pawcprj, only : pawcprj_type, pawcprj_alloc, pawcprj_free
51 : : use m_getghc, only : multithreaded_getghc
52 : : use m_gemm_nonlop_projectors , only : gemm_nonlop_use_gemm
53 : :
54 : : use m_xg
55 : : use m_xgTransposer
56 : :
57 : : #if defined(HAVE_GPU_MARKERS)
58 : : use m_nvtx_data
59 : : #endif
60 : :
61 : : #if defined(HAVE_GPU)
62 : : use m_gpu_toolbox
63 : : #endif
64 : :
65 : : #if defined(HAVE_YAKL)
66 : : use gator_mod
67 : : #endif
68 : :
69 : : use, intrinsic :: iso_c_binding, only: c_associated,c_loc,c_ptr,c_f_pointer,c_double,c_size_t
70 : :
71 : : use m_xmpi
72 : : use m_xomp
73 : : #ifdef HAVE_OPENMP
74 : : use omp_lib
75 : : #endif
76 : :
77 : : implicit none
78 : :
79 : : private
80 : :
81 : : integer, parameter :: l_tim_getghc=7
82 : : real(dp), parameter :: inv_sqrt2 = 1/sqrt2
83 : :
84 : : ! For use in getghc_gsc1
85 : : integer, save :: l_cpopt
86 : : integer, save :: l_icplx
87 : : integer, save :: l_nspinor
88 : : logical, save :: l_paw
89 : : integer, save :: l_prtvol
90 : : integer, save :: l_sij_opt
91 : : integer, save :: l_paral_kgb
92 : : integer, save :: l_useria
93 : : integer, save :: l_block_sliced
94 : : integer, save :: l_invovl_algo
95 : : real(dp), save :: l_invovl_precision
96 : :
97 : : type(mpi_type),pointer,save :: l_mpi_enreg
98 : : type(gs_hamiltonian_type),pointer,save :: l_gs_hamk
99 : :
100 : : integer, parameter :: DEBUG_ROWS = 5
101 : : integer, parameter :: DEBUG_COLUMNS = 5
102 : :
103 : : public :: chebfiwf2
104 : :
105 : : CONTAINS !========================================================================================
106 : : !!***
107 : :
108 : : !!****f* m_chebfiwf/chebfiwf2
109 : : !! NAME
110 : : !! chebfiwf2
111 : : !!
112 : : !! FUNCTION
113 : : !! This routine updates the whole wave functions set at a given k-point,
114 : : !! using the Chebfi method (2021 version using xG abstraction layer)
115 : : !!
116 : : !! INPUTS
117 : : !! dtset= input variables for this dataset
118 : : !! mpi_enreg= MPI-parallelisation information
119 : : !! nband= number of bands at this k point
120 : : !! npw= number of plane waves at this k point
121 : : !! nspinor= number of spinorial components of the wavefunctions
122 : : !! prtvol= control print volume and debugging
123 : : !!
124 : : !! OUTPUT
125 : : !! eig(nband)= eigenvalues (hartree) for all bands
126 : : !! enl_out(nband)= contribution of each band to the nl part of energy
127 : : !! resid(nband)= residuals for each band
128 : : !!
129 : : !! SIDE EFFECTS
130 : : !! cg(2,npw*nspinor*nband)= planewave coefficients of wavefunctions
131 : : !! gs_hamk <type(gs_hamiltonian_type)>=all data for the hamiltonian at k
132 : : !!
133 : : !! SOURCE
134 : :
135 : 11256 : subroutine chebfiwf2(cg,dtset,eig,occ,enl_out,gs_hamk,mpi_enreg,&
136 : 11256 : & nband,npw,nspinor,prtvol,resid)
137 : :
138 : : implicit none
139 : :
140 : : ! Arguments ------------------------------------
141 : : integer,intent(in) :: nband,npw,prtvol,nspinor
142 : : type(mpi_type),target,intent(in) :: mpi_enreg
143 : : real(dp),target,intent(inout) :: cg(2,npw*nspinor*nband)
144 : : real(dp),target,intent(out) :: resid(nband)
145 : : real(dp),intent(out) :: enl_out(nband)
146 : : real(dp),target,intent(out) :: eig(nband)
147 : : real(dp),target,intent(in) :: occ(nband)
148 : : type(dataset_type),intent(in) :: dtset
149 : : type(gs_hamiltonian_type),target,intent(inout) :: gs_hamk
150 : :
151 : : ! Local variables-------------------------------
152 : : ! scalars
153 : : integer, parameter :: tim_chebfiwf2 = 1750
154 : : integer, parameter :: tim_nonlop = 1753
155 : : integer :: iband,shift,space,blockdim,total_spacedim,ierr
156 : : integer :: me_g0,me_g0_fft
157 : : logical :: transfer_cg
158 : : integer(kind=c_size_t) :: localMem
159 : 11256 : type(chebfi_t) :: chebfi
160 : : type(xgBlock_t) :: xgx0,xgeigen,xgocc,xgresidu
161 : : ! arrays
162 : : real(dp) :: tsec(2)
163 : : integer(kind=c_size_t) :: chebfiMem(2)
164 : 11256 : real(dp), allocatable :: gvnlxc(:,:),occ_tmp(:)
165 : :
166 : : ! Parameters for nonlop call in NC
167 : : integer,parameter :: choice=1, paw_opt=0, signs=1
168 : : real(dp) :: gsc_dummy(1,1)
169 [ + + + + : 158224 : type(pawcprj_type) :: cprj_dum(gs_hamk%natom,1)
+ + + + +
- + - ]
170 : :
171 : : ! *********************************************************************
172 : :
173 : : !################ INITIALIZATION #####################################
174 : : !######################################################################
175 : :
176 : 11256 : call timab(tim_chebfiwf2,1,tsec)
177 : :
178 : : !Set module variables
179 : 11256 : l_paw = (gs_hamk%usepaw==1)
180 [ + + ]: 11256 : l_cpopt=-1;l_sij_opt=0;if (l_paw) l_sij_opt=1
181 : 11256 : l_nspinor = nspinor
182 : 11256 : l_prtvol = prtvol
183 : 11256 : l_mpi_enreg => mpi_enreg
184 : 11256 : l_gs_hamk => gs_hamk
185 : 11256 : l_paral_kgb = dtset%paral_kgb
186 : 11256 : l_invovl_algo = dtset%invovl_algo
187 : 11256 : l_block_sliced = dtset%invovl_blksliced
188 : 11256 : l_invovl_precision = dtset%invovl_precision
189 : :
190 : : !Variables
191 : 11256 : blockdim=l_mpi_enreg%nproc_band*l_mpi_enreg%bandpp
192 : : !for debug
193 : 11256 : l_useria=dtset%useria
194 : :
195 : : !Depends on istwfk
196 [ + + ]: 11256 : if ( gs_hamk%istwf_k > 1 ) then ! Real only
197 : : ! SPACE_CR mean that we have complex numbers but no re*im terms only re*re
198 : : ! and im*im so that a vector of complex is consider as a long vector of real
199 : : ! therefore the number of data is (2*npw*nspinor)*nband
200 : : ! This space is completely equivalent to SPACE_R but will correctly set and
201 : : ! get the array data into the xgBlock
202 : 3667 : space = SPACE_CR
203 : 3667 : l_icplx = 2
204 : : else ! complex
205 : 7589 : space = SPACE_C
206 : 7589 : l_icplx = 1
207 : : end if
208 : :
209 : 11256 : me_g0 = -1
210 : 11256 : me_g0_fft = -1
211 [ + + ]: 11256 : if (space==SPACE_CR) then
212 : 3667 : me_g0 = 0
213 : 3667 : me_g0_fft = 0
214 [ + + ]: 3667 : if (gs_hamk%istwf_k == 2) then
215 [ + + ]: 573 : if (l_mpi_enreg%me_g0 == 1) me_g0 = 1
216 [ + - ]: 573 : if (l_mpi_enreg%me_g0_fft == 1) me_g0_fft = 1
217 : : end if
218 : : end if
219 : :
220 : : !Memory info
221 [ - + ]: 11256 : if ( prtvol >= 3 ) then
222 [ # # ]: 0 : if (l_mpi_enreg%paral_kgb == 1) then
223 : 0 : total_spacedim = l_icplx*npw*nspinor
224 : 0 : call xmpi_sum(total_spacedim,l_mpi_enreg%comm_bandspinorfft,ierr)
225 : : else
226 : 0 : total_spacedim = 0
227 : : end if
228 : : chebfiMem = chebfi_memInfo(nband,l_icplx*npw*nspinor,space,l_mpi_enreg%paral_kgb, &
229 : 0 : & total_spacedim,l_mpi_enreg%bandpp) !blockdim
230 : 0 : localMem = (int(2,c_size_t)*npw*nspinor*nband+3*nband)*kind(1.d0) !blockdim
231 : 0 : write(std_out,'(1x,A,F10.6,1x,A)') "Each MPI process calling chebfi should need around ", &
232 [ # # ]: 0 : (localMem+sum(chebfiMem))/1e9,"GB of peak memory as follows :"
233 : 0 : write(std_out,'(4x,A,F10.6,1x,A)') "Permanent memory in chebfiwf : ",real(localMem)/1e9,"GB"
234 : 0 : write(std_out,'(4x,A,F10.6,1x,A)') "Permanent memory in m_chebfi : ",real(chebfiMem(1))/1e9,"GB"
235 : 0 : write(std_out,'(4x,A,F10.6,1x,A)') "Temporary memory in m_chebfi : ",real(chebfiMem(2))/1e9,"GB"
236 : : end if
237 : :
238 : 11256 : transfer_cg = .false.
239 : : #ifdef HAVE_OPENMP_OFFLOAD
240 : : !$OMP TARGET ENTER DATA MAP(to:eig,resid,occ) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
241 : : transfer_cg = .not. xomp_target_is_present(c_loc(cg))
242 : : !$OMP TARGET ENTER DATA MAP(to:cg) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP .and. transfer_cg)
243 : : #endif
244 : :
245 : : call xgBlock_map(xgx0,cg,space,npw*nspinor,nband,comm=l_mpi_enreg%comm_bandspinorfft,me_g0=me_g0,&
246 : 11256 : & gpu_option=dtset%gpu_option)
247 : :
248 : 11256 : call xgBlock_map_1d(xgeigen,eig,SPACE_R,nband,gpu_option=dtset%gpu_option)
249 : :
250 : 11256 : call xgBlock_map_1d(xgresidu,resid,SPACE_R,nband,gpu_option=dtset%gpu_option)
251 : :
252 : : ! Occupancies in chebyshev are used for convergence criteria only
253 [ + + + + : 11256 : if (dtset%nbdbuf==-101.and.nspinor==1.and.dtset%nsppol==1) then
+ - ]
254 [ - + + - : 2448 : ABI_MALLOC(occ_tmp,(nband))
+ - + - -
+ - + -
+ ]
255 [ + + ]: 7344 : occ_tmp(:) = half*occ(:)
256 : 816 : call xgBlock_map_1d(xgocc,occ_tmp,SPACE_R,nband,gpu_option=dtset%gpu_option)
257 : : else
258 : 10440 : call xgBlock_map_1d(xgocc,occ,SPACE_R,nband,gpu_option=dtset%gpu_option)
259 : : end if
260 : :
261 : 11256 : call timab(tim_chebfiwf2,2,tsec) ! IL TODO this should be deactivated
262 : :
263 : : ABI_NVTX_START_RANGE(NVTX_CHEBFI2_INIT)
264 : : call chebfi_init(chebfi,nband,npw*nspinor,dtset%tolwfr_diago,dtset%ecut, &
265 : : & dtset%paral_kgb,l_mpi_enreg%bandpp, &
266 : : & dtset%mdeg_filter, dtset%nbdbuf, space,1, &
267 : : & l_mpi_enreg%comm_bandspinorfft,me_g0,me_g0_fft,l_paw,&
268 : : & l_mpi_enreg%comm_spinorfft,l_mpi_enreg%comm_band,&
269 : : & dtset%chebfi_oracle,dtset%oracle_factor,dtset%oracle_min_occ,&
270 : : & l_gs_hamk%gpu_option,gpu_kokkos_nthrd=dtset%gpu_kokkos_nthrd,&
271 : 11256 : & gpu_thread_limit=dtset%gpu_thread_limit)
272 : : ABI_NVTX_END_RANGE()
273 : :
274 : : !################ RUUUUUUUN #####################################
275 : : !######################################################################
276 : :
277 : 11256 : call chebfi_run(chebfi,xgx0,getghc_gsc1,getBm1X,xgeigen,xgocc,xgresidu,nspinor)
278 : :
279 [ + + ]: 11256 : if (allocated(occ_tmp)) then
280 : 816 : ABI_FREE(occ_tmp)
281 : : end if
282 : :
283 [ + + ]: 11256 : if ( .not. l_paw ) then
284 : 1496 : call timab(tim_nonlop,1,tsec)
285 : : #ifdef FC_CRAY
286 : : ABI_MALLOC(gvnlxc,(1,1))
287 : : #else
288 [ - + - + ]: 1496 : ABI_MALLOC(gvnlxc,(0,0))
289 : : #endif
290 : : !end if
291 : :
292 : : ABI_NVTX_START_RANGE(NVTX_CHEBFI2_NONLOP)
293 : : !Call nonlop
294 [ + + ]: 1496 : if (l_paral_kgb==0) then
295 : :
296 : : call nonlop(choice,l_cpopt,cprj_dum,enl_out,l_gs_hamk,0,eig,mpi_enreg,nband,1,paw_opt,&
297 : 408 : & signs,gsc_dummy,l_tim_getghc,cg,gvnlxc)
298 : :
299 : : else
300 : : #ifdef HAVE_OPENMP_OFFLOAD
301 : : !$OMP TARGET UPDATE FROM(cg) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
302 : : #endif
303 [ + + ]: 2176 : do iband=1,nband/blockdim
304 : 1088 : shift = (iband-1)*blockdim*npw*nspinor
305 : : call prep_nonlop(choice,l_cpopt,cprj_dum, &
306 : : & enl_out((iband-1)*blockdim+1:iband*blockdim),l_gs_hamk,0,&
307 : : & eig((iband-1)*blockdim+1:iband*blockdim),blockdim,mpi_enreg,1,paw_opt,signs,&
308 : : & gsc_dummy,l_tim_getghc,cg(:,shift+1:shift+blockdim*npw*nspinor),gvnlxc(:,:),&
309 : 2176 : & already_transposed=.false.)
310 : : end do
311 : : end if
312 : : ABI_NVTX_END_RANGE()
313 [ - + ]: 1496 : ABI_FREE(gvnlxc)
314 : 1496 : call timab(tim_nonlop,2,tsec)
315 : : end if
316 : :
317 : : !Free chebfi
318 : 11256 : call chebfi_free(chebfi)
319 : :
320 : : #ifdef HAVE_OPENMP_OFFLOAD
321 : : !$OMP TARGET UPDATE FROM(eig,resid) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
322 : : !$OMP TARGET EXIT DATA MAP(delete:eig,resid,occ) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
323 : : !$OMP TARGET EXIT DATA MAP(from:cg) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP .and. transfer_cg)
324 : : #endif
325 : :
326 : 11256 : call timab(tim_chebfiwf2,2,tsec)
327 : :
328 : : DBG_EXIT("COLL")
329 : :
330 [ - + - + : 33768 : end subroutine chebfiwf2
- + - + -
+ ]
331 : : !!***
332 : :
333 : : !----------------------------------------------------------------------
334 : :
335 : : !!****f* m_chebfiwf/getghc_gsc1
336 : : !! NAME
337 : : !! getghc_gsc1
338 : : !!
339 : : !! FUNCTION
340 : : !! This routine computes H|C> and possibly S|C> for a given wave function C.
341 : : !! It acts as a driver for getghc, taken into account parallelism, multithreading, etc.
342 : : !!
343 : : !! SIDE EFFECTS
344 : : !! X <type(xgBlock_t)>= memory block containing |C>
345 : : !! AX <type(xgBlock_t)>= memory block containing H|C>
346 : : !! BX <type(xgBlock_t)>= memory block containing S|C>
347 : : !!
348 : : !! SOURCE
349 : :
350 : 92674 : subroutine getghc_gsc1(X,AX,BX)
351 : :
352 : : implicit none
353 : :
354 : : !Arguments ------------------------------------
355 : : type(xgBlock_t), intent(inout) :: X
356 : : type(xgBlock_t), intent(inout) :: AX
357 : : type(xgBlock_t), intent(inout) :: BX
358 : : integer :: blockdim
359 : : integer :: spacedim
360 [ + + + + : 1301916 : type(pawcprj_type) :: cprj_dum(l_gs_hamk%natom,1)
+ + + + +
- + - ]
361 : :
362 : : !Local variables-------------------------------
363 : : !scalars
364 : : real(dp) :: eval
365 : : !arrays
366 : 92674 : real(dp), pointer :: cg(:,:)
367 : 92674 : real(dp), pointer :: ghc(:,:)
368 : 92674 : real(dp), pointer :: gsc(:,:)
369 : : real(dp) :: gvnlxc(1,1)
370 : :
371 : : ! *********************************************************************
372 : :
373 : : ABI_NVTX_START_RANGE(NVTX_GETGHC)
374 : :
375 : 92674 : call xgBlock_getSize(X,spacedim,blockdim)
376 : 92674 : call xgBlock_check(X,AX)
377 : 92674 : call xgBlock_check(X,BX)
378 : :
379 : 92674 : call xgBlock_reverseMap(X,cg,rows=1,cols=spacedim*blockdim)
380 : 92674 : call xgBlock_reverseMap(AX,ghc,rows=1,cols=spacedim*blockdim)
381 : 92674 : call xgBlock_reverseMap(BX,gsc,rows=1,cols=spacedim*blockdim)
382 : :
383 : : call multithreaded_getghc(l_cpopt,cg,cprj_dum,ghc,gsc,&
384 : 92674 : l_gs_hamk,gvnlxc,eval,l_mpi_enreg,blockdim,l_prtvol,l_sij_opt,l_tim_getghc,0)
385 : :
386 : : #if defined(HAVE_GPU_CUDA) && defined(HAVE_YAKL)
387 : : call gpu_device_synchronize()
388 : : #endif
389 : :
390 [ + + ]: 92674 : if ( .not. l_paw ) call xgBlock_copy(X,BX)
391 : :
392 : : ABI_NVTX_END_RANGE()
393 : :
394 : 185348 : end subroutine getghc_gsc1
395 : : !!***
396 : :
397 : : !----------------------------------------------------------------------
398 : :
399 : : !!****f* m_chebfiwf/getBm1X
400 : : !! NAME
401 : : !! getBm1X
402 : : !!
403 : : !! FUNCTION
404 : : !! This routine computes S^-1|C> for a given wave function C.
405 : : !! It acts as a driver for apply_invovl.
406 : : !!
407 : : !! SIDE EFFECTS
408 : : !! X <type(xgBlock_t)>= memory block containing |C>
409 : : !! Bm1X <type(xgBlock_t)>= memory block containing S^-1|C>
410 : : !!
411 : : !! SOURCE
412 : :
413 : 69574 : subroutine getBm1X(X,Bm1X)
414 : :
415 : : implicit none
416 : :
417 : : !Arguments ------------------------------------
418 : : type(xgBlock_t), intent(inout) :: X
419 : : type(xgBlock_t), intent(inout) :: Bm1X
420 : :
421 : : !Local variables-------------------------------
422 : : !scalars
423 : : integer :: blockdim
424 : : integer :: spacedim
425 : : !arrays
426 : 69574 : real(dp), pointer :: ghc_filter(:,:)
427 : 69574 : real(dp), pointer :: gsm1hc_filter(:,:)
428 : 69574 : type(pawcprj_type), allocatable :: cwaveprj_next(:,:) !dummy
429 : :
430 : : ! *********************************************************************
431 : :
432 : 69574 : call xgBlock_getSize(X,spacedim,blockdim)
433 : :
434 [ + - ]: 69574 : if(l_paw) then
435 : :
436 : 69574 : call xgBlock_reverseMap(X,ghc_filter,rows=1,cols=spacedim*blockdim)
437 : 69574 : call xgBlock_reverseMap(Bm1X,gsm1hc_filter,rows=1,cols=spacedim*blockdim)
438 : :
439 : : !cwaveprj_next is dummy
440 [ + + ]: 69574 : if(gemm_nonlop_use_gemm) then
441 [ - + + + : 1152 : ABI_MALLOC(cwaveprj_next, (1,1))
+ + ]
442 : : else
443 [ - + + - : 1825960 : ABI_MALLOC(cwaveprj_next, (l_gs_hamk%natom,l_nspinor*blockdim))
- + + - +
- + - - +
- + + + +
+ ]
444 : 69286 : call pawcprj_alloc(cwaveprj_next,0,l_gs_hamk%dimcprj)
445 : : end if
446 : :
447 : : ABI_NVTX_START_RANGE(NVTX_INVOVL)
448 : : call apply_invovl(l_gs_hamk, ghc_filter(:,:), gsm1hc_filter(:,:), cwaveprj_next(:,:), &
449 [ - + - + ]: 69574 : spacedim/l_nspinor, blockdim, l_mpi_enreg, l_nspinor, l_block_sliced, l_invovl_precision, algo=l_invovl_algo)
450 : : ABI_NVTX_END_RANGE()
451 : :
452 : 69574 : call pawcprj_free(cwaveprj_next)
453 [ + - + + : 1117766 : ABI_FREE(cwaveprj_next)
- + - + -
+ ]
454 : :
455 : : else
456 : :
457 : 0 : call xgBlock_copy(X,Bm1X)
458 : :
459 : : end if
460 : :
461 [ - + - - : 139148 : end subroutine getBm1X
- - - - ]
462 : : !!***
463 : :
464 : : end module m_chebfiwf
465 : : !!***
|