Branch data Line data Source code
1 : : !!****f* ABINIT/m_lobpcgwf
2 : : !! NAME
3 : : !! m_lobpcgwf
4 : : !!
5 : : !! FUNCTION
6 : : !! This routine updates the whole wave functions at a given k-point,
7 : : !! using the lobpcg method
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) 1998-2026 ABINIT group (JB)
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_lobpcgwf
31 : :
32 : : use defs_basis
33 : : use m_abicore
34 : : use m_lobpcg
35 : : use m_xmpi
36 : : use m_errors
37 : : use m_time
38 : : use m_xomp
39 : : use m_fstrings
40 : : use m_xg
41 : : use m_xgTransposer
42 : : use m_lobpcg2
43 : : use m_dtset
44 : :
45 : : use defs_abitypes, only : mpi_type
46 : : use m_hamiltonian, only : gs_hamiltonian_type, K_H_K, KPRIME_H_KPRIME
47 : : use m_cgtools, only : cg_copy_spin
48 : : use m_pawcprj, only : pawcprj_type
49 : : use m_nonlop, only : nonlop
50 : : use m_prep_kgb, only : prep_getghc,prep_nonlop
51 : : use m_getghc, only : multithreaded_getghc
52 : :
53 : : #if defined(HAVE_GPU)
54 : : use m_gpu_toolbox
55 : : #endif
56 : :
57 : : #if defined(HAVE_GPU_MARKERS)
58 : : use m_nvtx_data
59 : : #endif
60 : :
61 : : use, intrinsic :: iso_c_binding
62 : :
63 : : implicit none
64 : : private
65 : :
66 : : integer, parameter :: l_tim_getghc=5
67 : : double precision, parameter :: inv_sqrt2 = 1/sqrt2
68 : :
69 : : ! For use in getghc_gsc1
70 : : integer,save :: l_cpopt
71 : : logical,save :: l_paw
72 : : integer,save :: l_prtvol
73 : : integer,save :: l_sij_opt
74 : : type(mpi_type),pointer,save :: l_mpi_enreg
75 : : type(gs_hamiltonian_type),pointer,save :: l_gs_hamk
76 : :
77 : : public :: lobpcgwf2
78 : :
79 : : contains
80 : :
81 [ - + ]: 24196 : subroutine lobpcgwf2(cg,dtset,eig,occ,enl_out,gs_hamk,isppol,ikpt,inonsc,istep,kinpw,kinpw_kp,mpi_enreg,&
82 : 24196 : & nband,npw,nspinor,prtvol,resid,nbdbuf)
83 : :
84 : :
85 : :
86 : : !Arguments ------------------------------------
87 : : integer,intent(in) :: nband,npw,prtvol,nspinor
88 : : integer,intent(in) :: isppol,ikpt,inonsc,istep,nbdbuf
89 : : type(gs_hamiltonian_type),target,intent(inout) :: gs_hamk
90 : : type(dataset_type) ,intent(in ) :: dtset
91 : : type(mpi_type) ,target,intent(in) :: mpi_enreg
92 : : real(dp) ,target,intent(inout) :: cg(2,nspinor*nband*npw)
93 : : real(dp) ,intent(in ) :: kinpw(npw),kinpw_kp(:)
94 : : real(dp) ,target,intent( out) :: resid(nband)
95 : : real(dp) ,intent( out) :: enl_out(nband)
96 : : real(dp) ,target,intent( out) :: eig(nband)
97 : : real(dp) ,target,intent(in ) :: occ(nband)
98 : :
99 : : !Local variables-------------------------------
100 : :
101 : : type(xgBlock_t) :: xgx0
102 : : type(xgBlock_t) :: xgeigen
103 : : type(xgBlock_t) :: xgresidu
104 : : type(xgBlock_t) :: xgocc
105 : : type(xgBlock_t) :: xg_precond
106 : 24196 : type(lobpcg_t) :: lobpcg
107 : :
108 : : integer :: space, blockdim,pcon_nspinor
109 : :
110 : : logical :: transfer_cg
111 : :
112 : : integer, parameter :: tim_lobpcgwf2 = 1640
113 : : integer, parameter :: tim_enl = 1657
114 : : double precision :: tsec(2)
115 : :
116 : : ! Important things for NC
117 : : integer,parameter :: choice=1, paw_opt=0, signs=1
118 [ + + + + : 193568 : type(pawcprj_type) :: cprj_dum(1,1)
+ + + - +
- ]
119 : : integer :: iblock, shift, me_g0, me_g0_fft
120 : : real(dp) :: gsc_dummy(0,0)
121 : 24196 : real(dp), allocatable :: gvnlxc(:,:)
122 : 24196 : real(dp), allocatable :: pcon(:),occ_tmp(:),cwavef_spin(:,:),enlout_spin(:)
123 : :
124 : : ! *********************************************************************
125 : :
126 : 24196 : call timab(tim_lobpcgwf2,1,tsec)
127 : :
128 : : ! Set module variables
129 : 24196 : l_paw = (gs_hamk%usepaw==1)
130 [ + + ]: 24196 : l_cpopt=-1;l_sij_opt=0;if (l_paw) l_sij_opt=1
131 : 24196 : l_prtvol = prtvol
132 : 24196 : l_mpi_enreg => mpi_enreg
133 : 24196 : l_gs_hamk => gs_hamk
134 : :
135 : : !Variables
136 : 24196 : blockdim=nband/dtset%nblock_lobpcg
137 [ - + ]: 24196 : if (blockdim/=mpi_enreg%nproc_band*mpi_enreg%bandpp) then ! without this check computation of enl_out can be wrong
138 : 0 : ABI_ERROR('blockdim is not consistent with nproc_band and bandpp')
139 : : end if
140 : :
141 : : !Depends on istwfk
142 [ + + ]: 24196 : if ( gs_hamk%istwf_k > 1 ) then ! Real only
143 : : ! SPACE_CR mean that we have complex numbers but no re*im terms only re*re
144 : : ! and im*im so that a vector of complex is consider as a long vector of real
145 : : ! therefore the number of data is (2*npw*nspinor)*nband
146 : : ! This space is completely equivalent to SPACE_R but will correctly set and
147 : : ! get the array data into the xgBlock
148 : 5036 : space = SPACE_CR
149 : : else ! complex
150 : 19160 : space = SPACE_C
151 : : end if
152 : :
153 : : !For preconditionning
154 [ + - ]: 24196 : if (gs_hamk%use_gbt == 0) then
155 [ - + + - : 72588 : ABI_MALLOC(pcon,(npw))
+ - + - -
+ - + -
+ ]
156 : 24196 : call build_pcon(pcon,kinpw,npw)
157 : 24196 : pcon_nspinor=nspinor
158 : : else
159 [ # # ]: 0 : ABI_CHECK(nspinor == 2,"GBT LOBPCG2 preconditioning requires two local spinor components")
160 [ # # ]: 0 : ABI_CHECK(size(kinpw_kp) >= npw,"GBT LOBPCG2 requires the k-prime kinetic table")
161 [ # # # # : 0 : ABI_MALLOC(pcon,(2*npw))
# # # # #
# # # # #
# # ]
162 : 0 : call build_pcon(pcon(1:npw),kinpw,npw)
163 : 0 : call build_pcon(pcon(npw+1:2*npw),kinpw_kp,npw)
164 : 0 : pcon_nspinor=1
165 : : end if
166 : :
167 : 24196 : transfer_cg = .false.
168 : : #ifdef HAVE_OPENMP_OFFLOAD
169 : : if(gs_hamk%gpu_option==ABI_GPU_OPENMP) then
170 : : !$OMP TARGET ENTER DATA MAP(to:eig,resid,occ,pcon)
171 : : transfer_cg = .not. xomp_target_is_present(c_loc(cg))
172 : : !$OMP TARGET ENTER DATA MAP(to:cg) IF(transfer_cg)
173 : : end if
174 : : #endif
175 : :
176 : : ! Local variables for lobpcg
177 : 24196 : me_g0 = -1
178 : 24196 : me_g0_fft = -1
179 [ + + ]: 24196 : if (space==SPACE_CR) then
180 : 5036 : me_g0 = 0
181 : 5036 : me_g0_fft = 0
182 [ + + ]: 5036 : if (gs_hamk%istwf_k == 2) then
183 [ + + ]: 2061 : if (mpi_enreg%me_g0 == 1) me_g0 = 1
184 [ + + ]: 2061 : if (mpi_enreg%me_g0_fft == 1) me_g0_fft = 1
185 : : end if
186 : : end if
187 : : call xgBlock_map(xgx0,cg,space,npw*nspinor,nband,comm=mpi_enreg%comm_bandspinorfft,me_g0=me_g0,&
188 : 24196 : & gpu_option=dtset%gpu_option)
189 : :
190 : 24196 : call xgBlock_map_1d(xg_precond,pcon,SPACE_R,npw*nspinor/pcon_nspinor,gpu_option=dtset%gpu_option)
191 : :
192 : 24196 : call xgBlock_map_1d(xgeigen,eig,SPACE_R,nband,gpu_option=dtset%gpu_option)
193 : :
194 : 24196 : call xgBlock_map_1d(xgresidu,resid,SPACE_R,nband,gpu_option=dtset%gpu_option)
195 : :
196 : : ! Occupancies in LOBPCG are used for convergence criteria only
197 [ + + + - : 24196 : if (dtset%nbdbuf==-101.and.nspinor==1.and.dtset%nsppol==1) then
- + ]
198 [ # # # # : 0 : ABI_MALLOC(occ_tmp,(nband))
# # # # #
# # # #
# ]
199 [ # # ]: 0 : occ_tmp(:) = half*occ(:)
200 : : #ifdef HAVE_OPENMP_OFFLOAD
201 : : !$OMP TARGET ENTER DATA MAP(to:occ_tmp) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
202 : : #endif
203 : 0 : call xgBlock_map_1d(xgocc,occ_tmp,SPACE_R,nband,gpu_option=dtset%gpu_option)
204 : : else
205 : 24196 : call xgBlock_map_1d(xgocc,occ,SPACE_R,nband,gpu_option=dtset%gpu_option)
206 : : end if
207 : :
208 : : call lobpcg_init(lobpcg,nband,npw*nspinor,blockdim,dtset%tolwfr_diago,dtset%nline,&
209 : : space,mpi_enreg%comm_bandspinorfft,dtset%paral_kgb,mpi_enreg%comm_spinorfft,mpi_enreg%comm_band,&
210 : 24196 : me_g0,me_g0_fft,gs_hamk%gpu_option,gpu_thread_limit=dtset%gpu_thread_limit)
211 : :
212 : : ! Run lobpcg
213 : : call lobpcg_run(lobpcg,xgx0,getghc_gsc1,xg_precond,xgeigen,xgocc,xgresidu,prtvol,nspinor,isppol,ikpt,inonsc,istep,&
214 : 24196 : nbdbuf,pcond_nspinor=pcon_nspinor)
215 : :
216 [ - + ]: 24196 : if (allocated(occ_tmp)) then
217 : : #ifdef HAVE_OPENMP_OFFLOAD
218 : : !$OMP TARGET EXIT DATA MAP(delete:occ_tmp) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
219 : : #endif
220 : 0 : ABI_FREE(occ_tmp)
221 : : end if
222 : : ! Free preconditionning since not needed anymore
223 : : #ifdef HAVE_OPENMP_OFFLOAD
224 : : !$OMP TARGET EXIT DATA MAP(delete:pcon) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
225 : : #endif
226 [ - + ]: 24196 : ABI_FREE(pcon)
227 : :
228 [ + + ]: 24196 : if ( .not. l_paw ) then
229 : 8820 : call timab(tim_enl,1,tsec)
230 : : #ifdef FC_CRAY
231 : : ABI_MALLOC(gvnlxc,(1,1))
232 : : #else
233 [ - + - + ]: 8820 : ABI_MALLOC(gvnlxc,(0,0))
234 : : #endif
235 : :
236 : : !Call nonlop
237 [ + + ]: 8820 : if (dtset%paral_kgb==0) then
238 : :
239 [ + - ]: 552 : if (gs_hamk%use_gbt == 0) then
240 : : call nonlop(choice,l_cpopt,cprj_dum,enl_out,l_gs_hamk,0,eig,mpi_enreg,nband,1,paw_opt,&
241 : 552 : & signs,gsc_dummy,l_tim_getghc,cg,gvnlxc)
242 : : else
243 : 0 : gs_hamk%nspinor=1
244 [ # # # # : 0 : ABI_MALLOC(cwavef_spin,(2,npw*nband))
# # # # #
# # # #
# ]
245 [ # # # # : 0 : ABI_MALLOC(enlout_spin,(nband))
# # # # #
# # # #
# ]
246 : :
247 : 0 : gs_hamk%ispin_gbt=1
248 : 0 : call cg_copy_spin(1,npw,nspinor,nband,cg,cwavef_spin)
249 : : call nonlop(choice,l_cpopt,cprj_dum,enl_out,l_gs_hamk,0,eig,mpi_enreg,nband,1,paw_opt,&
250 : 0 : & signs,gsc_dummy,l_tim_getghc,cwavef_spin,gvnlxc,select_k=K_H_K)
251 : :
252 : 0 : gs_hamk%ispin_gbt=2
253 : 0 : call cg_copy_spin(2,npw,nspinor,nband,cg,cwavef_spin)
254 : : call nonlop(choice,l_cpopt,cprj_dum,enlout_spin,l_gs_hamk,0,eig,mpi_enreg,nband,1,paw_opt,&
255 : 0 : & signs,gsc_dummy,l_tim_getghc,cwavef_spin,gvnlxc,select_k=KPRIME_H_KPRIME)
256 : :
257 [ # # ]: 0 : enl_out(:)=enl_out(:)+enlout_spin(:)
258 : 0 : gs_hamk%nspinor=2
259 [ # # ]: 0 : ABI_FREE(cwavef_spin)
260 [ # # ]: 0 : ABI_FREE(enlout_spin)
261 : : end if
262 : :
263 : : else
264 : : #ifdef HAVE_OPENMP_OFFLOAD
265 : : !$OMP TARGET UPDATE FROM(cg) IF(gs_hamk%gpu_option==ABI_GPU_OPENMP)
266 : : #endif
267 [ + - ]: 8268 : if (gs_hamk%use_gbt == 0) then
268 [ + + ]: 53964 : do iblock=1,nband/blockdim
269 : 45696 : shift = (iblock-1)*blockdim*npw*nspinor
270 : : call prep_nonlop(choice,l_cpopt,cprj_dum, &
271 : : & enl_out((iblock-1)*blockdim+1:iblock*blockdim),gs_hamk,0,&
272 : : & eig((iblock-1)*blockdim+1:iblock*blockdim),blockdim,mpi_enreg,1,paw_opt,signs,&
273 : : & gsc_dummy,l_tim_getghc,cg(:,shift+1:shift+blockdim*npw*nspinor),gvnlxc(:,:),&
274 : 53964 : & already_transposed=.false.)
275 : : end do
276 : : else
277 : 0 : gs_hamk%nspinor=1
278 [ # # # # : 0 : ABI_MALLOC(cwavef_spin,(2,npw*blockdim))
# # # # #
# # # #
# ]
279 [ # # # # : 0 : ABI_MALLOC(enlout_spin,(blockdim))
# # # # #
# # # #
# ]
280 [ # # ]: 0 : do iblock=1,nband/blockdim
281 : 0 : shift = (iblock-1)*blockdim*npw*nspinor
282 : :
283 : 0 : gs_hamk%ispin_gbt=1
284 : 0 : call cg_copy_spin(1,npw,nspinor,blockdim,cg(:,shift+1:shift+blockdim*npw*nspinor),cwavef_spin)
285 : : call prep_nonlop(choice,l_cpopt,cprj_dum,enl_out((iblock-1)*blockdim+1:iblock*blockdim),&
286 : : & gs_hamk,0,eig((iblock-1)*blockdim+1:iblock*blockdim),blockdim,mpi_enreg,1,paw_opt,signs,&
287 : 0 : & gsc_dummy,l_tim_getghc,cwavef_spin,gvnlxc(:,:),already_transposed=.false.,select_k=K_H_K)
288 : :
289 : 0 : gs_hamk%ispin_gbt=2
290 : 0 : call cg_copy_spin(2,npw,nspinor,blockdim,cg(:,shift+1:shift+blockdim*npw*nspinor),cwavef_spin)
291 : : call prep_nonlop(choice,l_cpopt,cprj_dum,enlout_spin,gs_hamk,0,&
292 : : & eig((iblock-1)*blockdim+1:iblock*blockdim),blockdim,mpi_enreg,1,paw_opt,signs,&
293 : 0 : & gsc_dummy,l_tim_getghc,cwavef_spin,gvnlxc(:,:),already_transposed=.false.,select_k=KPRIME_H_KPRIME)
294 : :
295 : : enl_out((iblock-1)*blockdim+1:iblock*blockdim)=&
296 [ # # ]: 0 : & enl_out((iblock-1)*blockdim+1:iblock*blockdim)+enlout_spin(:)
297 : : end do
298 : 0 : gs_hamk%nspinor=2
299 [ # # ]: 0 : ABI_FREE(cwavef_spin)
300 [ # # ]: 0 : ABI_FREE(enlout_spin)
301 : : end if
302 : : end if
303 [ - + ]: 8820 : ABI_FREE(gvnlxc)
304 : 8820 : call timab(tim_enl,2,tsec)
305 : : end if
306 : :
307 : : ! Free lobpcg
308 : 24196 : call lobpcg_free(lobpcg)
309 : :
310 : : #ifdef HAVE_OPENMP_OFFLOAD
311 : : if(gs_hamk%gpu_option==ABI_GPU_OPENMP) then
312 : : !$OMP TARGET EXIT DATA MAP(from:eig,resid,occ,pcon)
313 : : !$OMP TARGET EXIT DATA MAP(from:cg) IF(transfer_cg)
314 : : end if
315 : : #endif
316 : :
317 : 24196 : call timab(tim_lobpcgwf2,2,tsec)
318 : :
319 : : DBG_EXIT("COLL")
320 : :
321 [ - + - + : 48392 : end subroutine lobpcgwf2
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + ]
322 : :
323 : 709125 : subroutine getghc_gsc1(X,AX,BX)
324 : :
325 : : !Arguments ------------------------------------
326 : : type(xgBlock_t), intent(inout) :: X
327 : : type(xgBlock_t), intent(inout) :: AX
328 : : type(xgBlock_t), intent(inout) :: BX
329 : :
330 : : !Local variables-------------------------------
331 : : !scalars
332 : : integer :: blockdim
333 : : integer :: spacedim
334 : : real(dp) :: eval,dum
335 [ + + + + : 10500345 : type(pawcprj_type) :: cprj_dum(l_gs_hamk%natom,1)
+ + + + +
- + - ]
336 : : !arrays
337 : 709125 : real(dp), pointer :: cg(:,:)
338 : 709125 : real(dp), pointer :: ghc(:,:)
339 : 709125 : real(dp), pointer :: gsc(:,:)
340 : 709125 : real(dp), allocatable :: gvnlxc(:,:)
341 : :
342 : : ! *********************************************************************
343 : :
344 : : ABI_NVTX_START_RANGE(NVTX_GETGHC)
345 : :
346 : 709125 : call xgBlock_getSize(X,spacedim,blockdim)
347 : 709125 : call xgBlock_check(X,AX)
348 : 709125 : call xgBlock_check(X,BX)
349 : :
350 : 709125 : call xgBlock_reverseMap(X,cg,rows=1,cols=spacedim*blockdim)
351 : 709125 : call xgBlock_reverseMap(AX,ghc,rows=1,cols=spacedim*blockdim)
352 : 709125 : call xgBlock_reverseMap(BX,gsc,rows=1,cols=spacedim*blockdim)
353 : :
354 : : #ifdef FC_CRAY
355 : : ABI_MALLOC(gvnlxc,(1,1))
356 : : #else
357 [ - + ]: 709125 : ABI_MALLOC(gvnlxc,(0,0))
358 : : #endif
359 : :
360 [ + + + + ]: 709125 : if (l_mpi_enreg%nproc_fft==1.or.l_gs_hamk%istwf_k==1) then
361 : : call multithreaded_getghc(l_cpopt,cg,cprj_dum,ghc,gsc,&
362 : 707469 : l_gs_hamk,gvnlxc,eval,l_mpi_enreg,blockdim,l_prtvol,l_sij_opt,l_tim_getghc,0)
363 [ + - ]: 1656 : else if (l_gs_hamk%istwf_k==2) then ! nproc_fft>1 and istwfk==2
364 : : call prep_getghc(cg(:,1:blockdim*spacedim),l_gs_hamk,gvnlxc,ghc,gsc(:,1:blockdim*spacedim),dum,blockdim,&
365 : 1656 : l_mpi_enreg,l_prtvol,l_sij_opt,l_cpopt,cprj_dum,already_transposed=.true.)
366 : : else ! nproc_fft>1 and istwfk>2
367 : 0 : ABI_ERROR('getghc in lobpcg not implemented for npfft>1 and istwfk>2')
368 : : end if
369 : :
370 [ - + ]: 709125 : ABI_FREE(gvnlxc)
371 : :
372 : : #if defined(HAVE_GPU_CUDA) && defined(HAVE_YAKL)
373 : : call gpu_device_synchronize()
374 : : #endif
375 : :
376 [ + + ]: 709125 : if ( .not. l_paw ) call xgBlock_copy(X,BX)
377 : :
378 : : ABI_NVTX_END_RANGE()
379 : :
380 : 2127375 : end subroutine getghc_gsc1
381 : :
382 [ - + - + ]: 24196 : subroutine build_pcon(pcon,kinpw,npw)
383 : :
384 : : integer,intent(in) :: npw
385 : : real(dp),intent(in) :: kinpw(:)
386 : : real(dp),intent(out) :: pcon(:)
387 : :
388 : : integer :: ipw
389 : :
390 : : !$omp parallel do schedule(static), shared(pcon,kinpw)
391 [ + + ]: 5117830 : do ipw=1,npw
392 [ + + ]: 5117830 : if(kinpw(ipw)>huge(0.0_dp)*1.d-11) then
393 : 194498 : pcon(ipw)=0.d0
394 : : else
395 : : pcon(ipw) = (27+kinpw(ipw)*(18+kinpw(ipw)*(12+8*kinpw(ipw)))) &
396 : 4899136 : & / (27+kinpw(ipw)*(18+kinpw(ipw)*(12+8*kinpw(ipw))) + 16*kinpw(ipw)**4)
397 : : end if
398 : : end do
399 : :
400 : 24196 : end subroutine build_pcon
401 : :
402 : : end module m_lobpcgwf
403 : : !!***
|
