Line data Source code
1 : !!****f* ABINIT/m_chebfiwf_cprj
2 : !! NAME
3 : !! m_chebfiwf_cprj
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) 2023-2026 ABINIT group (LB)
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 : module m_chebfiwf_cprj
28 :
29 : use, intrinsic :: iso_c_binding
30 : use defs_abitypes
31 : use defs_basis
32 : use m_abicore
33 : use m_errors
34 : use m_fstrings
35 : use m_time
36 : use m_xg
37 : use m_xg_nonlop
38 : use m_chebfi2_cprj
39 : use m_chebfi
40 : use m_invovl
41 :
42 : use m_dtset, only : dataset_type
43 : use m_hamiltonian, only : gs_hamiltonian_type
44 : use m_pawcprj, only : pawcprj_type
45 : use m_getghc, only : multithreaded_getghc
46 :
47 : use m_xg
48 : use m_xgTransposer
49 :
50 : use iso_c_binding, only: c_associated,c_loc,c_ptr,c_f_pointer
51 :
52 : use m_xmpi
53 : use m_xomp
54 : #ifdef HAVE_OPENMP
55 : use omp_lib
56 : #endif
57 :
58 : implicit none
59 :
60 : private
61 :
62 : integer, parameter :: l_tim_getghc=5
63 :
64 : ! For use in getghc_gsc1
65 : integer, save :: l_prtvol
66 : type(mpi_type),pointer,save :: l_mpi_enreg
67 : type(gs_hamiltonian_type),pointer,save :: l_gs_hamk
68 :
69 : public :: chebfiwf2_cprj
70 :
71 : contains
72 : !!***
73 :
74 : !!****f* m_chebfiwf/chebfiwf2_cprj
75 : !! NAME
76 : !! chebfiwf2_cprj
77 : !!
78 : !! FUNCTION
79 : !! This routine updates the whole wave functions set at a given k-point,
80 : !! using the Chebfi method (2021 version using xG abstraction layer)
81 : !!
82 : !! INPUTS
83 : !! dtset= input variables for this dataset
84 : !! mpi_enreg= MPI-parallelisation information
85 : !! nband= number of bands at this k point
86 : !! npw= number of plane waves at this k point
87 : !! nspinor= number of spinorial components of the wavefunctions
88 : !! prtvol= control print volume and debugging
89 : !!
90 : !! OUTPUT
91 : !! eig(nband)= eigenvalues (hartree) for all bands
92 : !! enl_out(nband)= contribution of each band to the nl part of energy
93 : !! resid(nband)= residuals for each band
94 : !!
95 : !! SIDE EFFECTS
96 : !! cg(2,npw*nspinor*nband)= planewave coefficients of wavefunctions
97 : !! gs_hamk <type(gs_hamiltonian_type)>=all data for the hamiltonian at k
98 : !!
99 : !! SOURCE
100 :
101 15368 : subroutine chebfiwf2_cprj(cg,dtset,eig,occ,enl_out,gs_hamk,mpi_enreg,&
102 15368 : & nband,npw,nspinor,prtvol,resid,xg_nonlop)
103 :
104 : !Arguments ------------------------------------
105 : integer,intent(in) :: nband,npw,prtvol,nspinor
106 : type(gs_hamiltonian_type),target,intent(inout) :: gs_hamk
107 : type(dataset_type) ,intent(in ) :: dtset
108 : type(mpi_type) ,target,intent(in) :: mpi_enreg
109 : real(dp) ,target,intent(inout) :: cg(2,nspinor*nband*npw)
110 : real(dp) ,target,intent( out) :: resid(nband)
111 : real(dp) ,intent( out) :: enl_out(nband)
112 : real(dp) ,target,intent( out) :: eig(nband)
113 : real(dp) ,target,intent(in ) :: occ(nband)
114 : type(xg_nonlop_t) , intent(in ) :: xg_nonlop
115 :
116 : !Local variables-------------------------------
117 :
118 : type(xgBlock_t) :: xgx0
119 : type(xg_t) :: cprj_xgx0
120 : type(xgBlock_t) :: xgeigen
121 : type(xgBlock_t) :: xgocc
122 : type(xgBlock_t) :: xgresidu
123 : type(xgBlock_t) :: xgenl
124 : type(xgBlock_t) :: xg_kin
125 15368 : type(chebfi_t) :: chebfi
126 :
127 : logical :: paw
128 :
129 : integer :: space, space_cprj, cprjdim, nband_cprj
130 : integer :: me_g0,me_g0_fft
131 :
132 : integer, parameter :: tim_chebfiwf2 = 2060
133 : double precision :: tsec(2)
134 :
135 15368 : real(dp), allocatable :: kin(:),occ_tmp(:)
136 :
137 : ! *********************************************************************
138 :
139 15368 : call timab(tim_chebfiwf2,1,tsec)
140 :
141 15368 : paw = gs_hamk%usepaw==1
142 :
143 : ! Set module variables
144 15368 : l_prtvol = prtvol
145 15368 : l_mpi_enreg => mpi_enreg
146 15368 : l_gs_hamk => gs_hamk
147 :
148 15368 : cprjdim = xg_nonlop%cprjdim
149 15368 : nband_cprj=nband/mpi_enreg%nproc_band
150 :
151 : !Depends on istwfk
152 15368 : if ( gs_hamk%istwf_k > 1 ) then ! Real only
153 : ! SPACE_CR mean that we have complex numbers but no re*im terms only re*re
154 : ! and im*im so that a vector of complex is consider as a long vector of real
155 : ! therefore the number of data is (2*npw*nspinor)*nband
156 : ! This space is completely equivalent to SPACE_R but will correctly set and
157 : ! get the array data into the xgBlock
158 6256 : space = SPACE_CR
159 : else ! complex
160 9112 : space = SPACE_C
161 : end if
162 15368 : space_cprj = xg_nonlop%space_cprj
163 :
164 : !For kinetic part of the Hamiltonian
165 46104 : ABI_MALLOC(kin,(gs_hamk%npw_fft_k))
166 15368 : call build_kin(kin,gs_hamk%kinpw_k,gs_hamk%npw_fft_k)
167 15368 : call xgBlock_map_1d(xg_kin,kin,SPACE_R,gs_hamk%npw_fft_k)
168 :
169 : ! Local variables for chebfi
170 15368 : me_g0 = -1
171 15368 : me_g0_fft = -1
172 15368 : if (space==SPACE_CR) then
173 6256 : me_g0 = 0
174 6256 : me_g0_fft = 0
175 6256 : if (gs_hamk%istwf_k == 2) then
176 782 : if (l_mpi_enreg%me_g0 == 1) me_g0 = 1
177 782 : if (l_mpi_enreg%me_g0_fft == 1) me_g0_fft = 1
178 : end if
179 : end if
180 15368 : call xgBlock_map(xgx0,cg,space,npw*nspinor,nband,l_mpi_enreg%comm_band,me_g0=me_g0)
181 :
182 15368 : call xgBlock_map_1d(xgeigen,eig,SPACE_R,nband)
183 :
184 15368 : call xgBlock_map_1d(xgresidu,resid,SPACE_R,nband)
185 :
186 15368 : call xgBlock_map_1d(xgenl,enl_out,SPACE_R,nband)
187 :
188 : ! Occupancies in chebyshev are used for convergence criteria only
189 15368 : if (dtset%nbdbuf==-101.and.nspinor==1.and.dtset%nsppol==1) then
190 2448 : ABI_MALLOC(occ_tmp,(nband))
191 7344 : occ_tmp(:) = half*occ(:)
192 816 : call xgBlock_map_1d(xgocc,occ_tmp,SPACE_R,nband,gpu_option=dtset%gpu_option)
193 : else
194 14552 : call xgBlock_map_1d(xgocc,occ,SPACE_R,nband,gpu_option=dtset%gpu_option)
195 : end if
196 :
197 : !call xg_cprj_copy(cprj_cwavef_bands,cprj_contiguous,space_cprj,nband_cprj,cprj_xgx0,&
198 : ! & xg_nonlop,l_mpi_enreg%comm_band,CPRJ_ALLOC)
199 15368 : call xg_init(cprj_xgx0,space_cprj,xg_nonlop%cprjdim,nband_cprj*nspinor,comm=l_mpi_enreg%comm_band)
200 :
201 : call chebfi_init(chebfi,nband,npw*nspinor,cprjdim,dtset%tolwfr_diago,dtset%ecut, &
202 : & mpi_enreg%bandpp, dtset%nline, dtset%nbdbuf, space,space_cprj,1, &
203 : & l_mpi_enreg%comm_band,me_g0,paw,&
204 : & dtset%chebfi_oracle,dtset%oracle_factor,dtset%oracle_min_occ,&
205 15368 : & xg_nonlop,me_g0_fft)
206 :
207 :
208 : ! Run chebfi
209 15368 : call chebfi_run_cprj(chebfi,xgx0,cprj_xgx0%self,xg_getghc,xg_kin,xgeigen,xgocc,xgresidu,xgenl,nspinor)
210 :
211 15368 : if (allocated(occ_tmp)) then
212 816 : ABI_FREE(occ_tmp)
213 : end if
214 15368 : ABI_FREE(kin)
215 :
216 : ! call xg_cprj_copy(cprj_cwavef_bands,cprj_contiguous,space_cprj,nband_cprj,cprj_xgx0,&
217 : ! & xg_nonlop,l_mpi_enreg%comm_band,CPRJ_FREE)
218 15368 : call xg_free(cprj_xgx0)
219 :
220 : ! Free chebfi
221 15368 : call chebfi_free(chebfi)
222 :
223 15368 : call timab(tim_chebfiwf2,2,tsec)
224 :
225 : DBG_EXIT("COLL")
226 :
227 30736 : end subroutine chebfiwf2_cprj
228 : !!***
229 :
230 : !!****f* m_chebfi/xg_getghc
231 : !! NAME
232 : !! xg_getghc
233 : !!
234 : !! FUNCTION
235 : !! This routine computes H|C> and possibly S|C> for a given wave function C.
236 : !! It acts as a driver for getghc, taken into account parallelism, multithreading, etc.
237 : !!
238 : !! SIDE EFFECTS
239 : !! X <type(xgBlock_t)>= memory block containing |C>
240 : !! AX <type(xgBlock_t)>= memory block containing H|C>
241 : !!
242 : !! SOURCE
243 : !
244 129554 : subroutine xg_getghc(X,AX)
245 :
246 : !Arguments ------------------------------------
247 : type(xgBlock_t), intent(inout) :: X
248 : type(xgBlock_t), intent(inout) :: AX
249 :
250 : !Local variables-------------------------------
251 : !scalars
252 : integer :: blockdim
253 : integer :: spacedim
254 : integer,parameter :: sij_opt=0,cpopt=-1,type_calc=1 ! Compute local part only
255 : ! integer :: iatom,iband,ispinor,cprj_index,cprj_rows,cprj_cols,ncpgr,nlmn
256 : real(dp) :: eval
257 1818236 : type(pawcprj_type) :: cprj_dum(l_gs_hamk%natom,1)
258 : !arrays
259 129554 : real(dp), pointer :: cg(:,:)
260 129554 : real(dp), pointer :: ghc(:,:)
261 : real(dp) :: gsc(1,1),gvnlxc(1,1)
262 :
263 : ! *********************************************************************
264 :
265 129554 : call xgBlock_getSize(X,spacedim,blockdim)
266 129554 : call xgBlock_check(X,AX)
267 :
268 129554 : call xgBlock_reverseMap(X,cg,rows=1,cols=spacedim*blockdim)
269 129554 : call xgBlock_reverseMap(AX,ghc,rows=1,cols=spacedim*blockdim)
270 :
271 : ! Apply only local part of the Hamiltonian
272 : call multithreaded_getghc(cpopt,cg,cprj_dum,ghc,gsc,&
273 129554 : l_gs_hamk,gvnlxc,eval,l_mpi_enreg,blockdim,l_prtvol,sij_opt,l_tim_getghc,type_calc)
274 :
275 259108 : end subroutine xg_getghc
276 : !!***
277 :
278 15368 : subroutine build_kin(kin,kinpw,npw)
279 :
280 : integer,intent(in) :: npw
281 : real(dp),intent(in) :: kinpw(:)
282 : real(dp),intent(out) :: kin(:)
283 :
284 : integer :: ipw
285 :
286 : !$omp parallel do schedule(static), shared(kin,kinpw)
287 1692923 : do ipw=1,npw
288 1692923 : if(kinpw(ipw)>huge(0.0_dp)*1.d-11) then
289 43860 : kin(ipw)=0.d0
290 : else
291 1633695 : kin(ipw) = kinpw(ipw)
292 : end if
293 : end do
294 :
295 15368 : end subroutine build_kin
296 :
297 : end module m_chebfiwf_cprj
298 : !!***
|
