Branch data Line data Source code
1 : : !!****m* ABINIT/m_gwls_wf
2 : : !! NAME
3 : : !! m_gwls_wf
4 : : !!
5 : : !! FUNCTION
6 : : !! .
7 : : !!
8 : : !! COPYRIGHT
9 : : !! Copyright (C) 2009-2026 ABINIT group (JLJ, BR, MC)
10 : : !! This file is distributed under the terms of the
11 : : !! GNU General Public License, see ~abinit/COPYING
12 : : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : : !!
14 : : !! SOURCE
15 : :
16 : : #if defined HAVE_CONFIG_H
17 : : #include "config.h"
18 : : #endif
19 : :
20 : : #include "abi_common.h"
21 : :
22 : :
23 : : !---------------------------------------------------------------------
24 : : ! Modules to handle low-level Abinit entities, like the Hamiltonian
25 : : ! and wavefunctions.
26 : : !---------------------------------------------------------------------
27 : :
28 : : module m_gwls_wf
29 : :
30 : : ! local modules
31 : : use m_gwls_utility
32 : :
33 : : ! abinit modules
34 : : use defs_basis
35 : : use m_abicore
36 : : use m_cgtools
37 : : use m_xmpi
38 : :
39 : : implicit none
40 : : save
41 : : private
42 : : !!***
43 : :
44 : : real(dp) :: dv=0.
45 : : integer :: nx=0, ny=0, nz=0, nk=0, nb=0, ng=0, i=0, cbf=1, cf=1, cb=1
46 : : !integer :: j=0, k=0
47 : : !!***
48 : :
49 : : public :: set_wf, norm_k, norm_kc, scprod_kc, scprod_k
50 : : ! public :: contribution, contribution_bloc
51 : : !!***
52 : : contains
53 : :
54 : : !!****f* m_hamiltonian/set_wf
55 : : !! NAME
56 : : !! set_wf
57 : : !!
58 : : !! FUNCTION
59 : : !! .
60 : : !!
61 : : !! INPUTS
62 : : !!
63 : : !! OUTPUT
64 : : !!
65 : : !! SOURCE
66 : :
67 : 8 : subroutine set_wf(dv2,nx2,ny2,nz2,nk2,nb2,ng2,cbf2,cf2,cb2)
68 : :
69 : : real(dp), intent(in) :: dv2
70 : : integer, intent(in) :: nx2, ny2, nz2, nk2, nb2, ng2, cbf2, cf2, cb2
71 : : ! *************************************************************************
72 : 8 : dv=dv2
73 : 8 : nx=nx2
74 : 8 : ny=ny2
75 : 8 : nz=nz2
76 : 8 : nk=nk2
77 : 8 : nb=nb2
78 : 8 : ng=ng2
79 : 8 : cbf=cbf2
80 : 8 : cf=cf2
81 : 8 : cb=cb2
82 : 8 : end subroutine set_wf
83 : : !!***
84 : :
85 : : !!****f* m_hamiltonian/norm_k
86 : : !! NAME
87 : : !! norm_k
88 : : !!
89 : : !! FUNCTION
90 : : !! .
91 : : !!
92 : : !! INPUTS
93 : : !!
94 : : !! OUTPUT
95 : : !!
96 : : !!
97 : : !!
98 : : !! SOURCE
99 : :
100 : 480 : real(dp) function norm_k(v)
101 : :
102 : : real(dp), intent(in) :: v(2,nk)
103 : : ! *************************************************************************
104 : 240 : norm_k = 0.0_dp
105 : :
106 [ + + ]: 38790 : do i=1,nk
107 : 38790 : norm_k = norm_k + v(1,i)**2 + v(2,i)**2
108 : : end do
109 : 240 : call xmpi_sum(norm_k,cbf,i) ! sum on all processors
110 : :
111 : 240 : norm_k = dsqrt(norm_k)
112 : :
113 : 240 : end function norm_k
114 : : !!***
115 : :
116 : : !!****f* m_hamiltonian/norm_kc
117 : : !! NAME
118 : : !! norm_kc
119 : : !!
120 : : !! FUNCTION
121 : : !! .
122 : : !!
123 : : !! INPUTS
124 : : !!
125 : : !! OUTPUT
126 : : !!
127 : : !!
128 : : !!
129 : : !! SOURCE
130 : :
131 : 22440 : real(dp) function norm_kc(v)
132 : :
133 : : complex(dp), intent(in) :: v(nk)
134 : : ! *************************************************************************
135 : 11220 : norm_kc = zero
136 : :
137 [ + + ]: 2894760 : do i=1,nk
138 : 2894760 : norm_kc = norm_kc + dble(v(i))**2+dimag(v(i))**2
139 : : end do
140 : 11220 : call xmpi_sum(norm_kc,cbf,i) ! sum on all processors
141 : 11220 : norm_kc = dsqrt(norm_kc)
142 : :
143 : 11220 : end function norm_kc
144 : : !!***
145 : :
146 : : !!****f* m_hamiltonian/scprod_kc
147 : : !! NAME
148 : : !! scprod_kc
149 : : !!
150 : : !! FUNCTION
151 : : !! .
152 : : !!
153 : : !! INPUTS
154 : : !!
155 : : !! OUTPUT
156 : : !!
157 : : !!
158 : : !!
159 : : !! SOURCE
160 : :
161 : 14448 : complex(dp) function scprod_kc(v1,v2)
162 : :
163 : : complex(dp), intent(in) :: v1(nk), v2(nk)
164 : : ! *************************************************************************
165 : 7224 : scprod_kc = zero
166 [ + + ]: 1863792 : do i=1,nk
167 : 1863792 : scprod_kc = scprod_kc + conjg(v1(i))*v2(i)
168 : : end do
169 : 7224 : call xmpi_sum(scprod_kc,cbf,i) ! sum on all processors
170 : :
171 : : !scprod_kc = sum(conjg(v1)*v2) !Eliminated to avoid functions that return vectors embeded in another function/subroutine call.
172 : :
173 : 7224 : end function scprod_kc
174 : : !!***
175 : :
176 : : !!****f* m_hamiltonian/contribution
177 : : !! NAME
178 : : !! contribution
179 : : !!
180 : : !! FUNCTION
181 : : !! .
182 : : !!
183 : : !! INPUTS
184 : : !!
185 : : !! OUTPUT
186 : : !!
187 : : !!
188 : : !!
189 : : !! SOURCE
190 : :
191 : : ! real(dp) function contribution(alpha,beta,k,norm_svne)
192 : : !
193 : : ! real(dp), intent(in) :: alpha(k), beta(k-1), norm_svne
194 : : ! integer, intent(in) :: k
195 : : !
196 : : ! integer :: i
197 : : ! real(dp), allocatable :: eq_lin(:,:)
198 : : !! *************************************************************************
199 : : ! ABI_MALLOC(eq_lin,(4,k))
200 : : !
201 : : ! !Copy the input data into 4 vectors that will be overwritten by dgtsv()
202 : : ! eq_lin = zero
203 : : ! eq_lin(1,1:k-1) = beta !sub-diagonal (only the 1:kmax-1 elements are used)
204 : : ! eq_lin(2,:) = alpha !diagonal
205 : : ! eq_lin(3,1:k-1) = beta !supra-diagonal (only the 1:kmax-1 elements are used)
206 : : ! eq_lin(4,1) = 1.0 !the RHS vector to the linear equation (here, |1,0,0,...>)
207 : : !
208 : : ! !DGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
209 : : ! !dgtsv(matrix size, # of column of RHS, sub-diagonal elements, diagonal elements, super-diagonal elements,
210 : : ! ! RHS of equation (solution of equation at end of routine), size of RHS vector,
211 : : ! ! error message integer (0: success, -i: illegal ith argument, i: ith factorization failed))
212 : : ! call dgtsv(k,1,eq_lin(1,1:k-1),eq_lin(2,:),eq_lin(3,1:k-1),eq_lin(4,:),k,i)
213 : : ! !Do the scalar product between the <qr_k|sv|ne>=(1,0,0,...) vector and the solution x to T*x=(1,0,0,...) to obtain contribution
214 : : ! !to SEX, at order k, from orbital n.
215 : : ! !We now replace the screened coulomb interaction by the coulomb hole...
216 : : ! contribution = (eq_lin(4,1)-1.0)*norm_svne**2
217 : : ! ABI_FREE(eq_lin)
218 : : ! end function contribution
219 : : !!***
220 : :
221 : : !!****f* m_hamiltonian/contribution_bloc
222 : : !! NAME
223 : : !! contribution_bloc
224 : : !!
225 : : !! FUNCTION
226 : : !! .
227 : : !!
228 : : !! INPUTS
229 : : !!
230 : : !! OUTPUT
231 : : !!
232 : : !!
233 : : !!
234 : : !! SOURCE
235 : :
236 : : ! function contribution_bloc(alpha,beta,kmax,norm_svne,nseeds)
237 : : !
238 : : ! real(dp) :: contribution_bloc(2)
239 : : ! integer, intent(in) :: kmax, nseeds
240 : : ! real(dp), intent(in) :: alpha(2,nseeds,nseeds,kmax), beta(2,nseeds,nseeds,kmax-1), norm_svne
241 : : !
242 : : ! integer :: i=0, j=0, k=0
243 : : ! integer, allocatable :: ipiv(:)
244 : : ! complex(dp), allocatable :: a(:,:),b(:,:) !For the (non-banded) solver of AX=B
245 : : !! *************************************************************************
246 : : ! !write(std_out,*) "Allocating..."
247 : : ! ABI_MALLOC(a,(kmax*nseeds,kmax*nseeds))
248 : : ! ABI_MALLOC(b,(kmax*nseeds,1))
249 : : ! ABI_MALLOC(ipiv,(kmax*nseeds))
250 : : ! !write(std_out,*) "Zeroing..."
251 : : ! a=zero
252 : : ! b=zero
253 : : ! ipiv=zero
254 : : !
255 : : ! !Copy the input data into 4 vectors that will be overwritten by dgtsv()
256 : : ! !do j=1,k*nseeds
257 : : ! ! do i=max(1,j-nseeds),min(k*nseeds,j+nseeds)
258 : : ! ! ab(2*nseeds+1+i-j,j) = cmplx(alpha
259 : : ! ! end do
260 : : ! !end do
261 : : ! !write(std_out,*) "Copying alpha..."
262 : : ! do k=1,kmax
263 : : ! do j=1,nseeds
264 : : ! do i=1,nseeds
265 : : ! a((k-1)*nseeds+i,(k-1)*nseeds+j) = cmplx(alpha(1,i,j,k),alpha(2,i,j,k),dp)
266 : : ! end do
267 : : ! end do
268 : : ! end do
269 : : ! !write(std_out,*) "Copying beta..."
270 : : ! do k=1,kmax-1
271 : : ! do j=1,nseeds
272 : : ! do i=1,j
273 : : ! a(k*nseeds+i,(k-1)*nseeds+j) = cmplx(beta(1,i,j,k),beta(2,i,j,k),dp)
274 : : ! a((k-1)*nseeds+i,k*nseeds+j) = cmplx(beta(1,j,i,k),-beta(2,j,i,k),dp)
275 : : ! end do
276 : : ! end do
277 : : ! end do
278 : : ! !write(std_out,*) "Setting RHS..."
279 : : ! b(1,1) = (1.0,0.0)
280 : : !
281 : : ! !write(std_out,*) "Solving..."
282 : : ! call zgesv(kmax*nseeds,1,a,kmax*nseeds,ipiv,b,kmax*nseeds,i)
283 : : ! !Do the scalar product between the <qr_k|sv|ne>=(1,0,0,...) vector and the solution x to T*x=(1,0,0,...) to obtain contribution
284 : : ! !to SEX, at order k, from orbital n
285 : : ! !write(std_out,*) "Obtaining contribution..."
286 : : ! contribution_bloc(1) = real(b(1,1))*norm_svne**2
287 : : ! contribution_bloc(2) = aimag(b(1,1))*norm_svne**2
288 : : ! !write(std_out,*) "Deallocating..."
289 : : ! ABI_FREE(a)
290 : : ! ABI_FREE(b)
291 : : ! ABI_FREE(ipiv)
292 : : ! end function contribution_bloc
293 : : !!***
294 : :
295 : : !!****f* m_hamiltonian/scprod_k
296 : : !! NAME
297 : : !! scprod_k
298 : : !!
299 : : !! FUNCTION
300 : : !! .
301 : : !!
302 : : !! INPUTS
303 : : !!
304 : : !! OUTPUT
305 : : !!
306 : : !!
307 : : !!
308 : : !! SOURCE
309 : :
310 [ - + ]: 80 : function scprod_k(v1,v2)
311 : : !--------------------------------------------------------------------------------
312 : : ! This function computes the inner product of two "vectors" (typically
313 : : ! wavefunctions), < v1 | v2 >.
314 : : !--------------------------------------------------------------------------------
315 : : real(dp) :: scprod_k(2)
316 : : real(dp), intent(in) :: v1(2,nk), v2(2,nk)
317 : : ! *************************************************************************
318 [ + + ]: 240 : scprod_k = zero
319 : :
320 : : ! The mitaine way
321 : : ! do i=1,nk
322 : : ! scprod_k(1) = scprod_k(1) + v1(1,i)*v2(1,i) + v1(2,i)*v2(2,i)
323 : : ! scprod_k(2) = scprod_k(2) + v1(1,i)*v2(2,i) - v1(2,i)*v2(1,i)
324 : : ! end do
325 : :
326 : : ! The ABINIT way
327 : 80 : scprod_k = cg_zdotc(nk,v1,v2)
328 : :
329 : : ! Collect from every processor
330 : 80 : call xmpi_sum(scprod_k,cbf,i) ! sum on all processors
331 : :
332 : : end function scprod_k
333 : : !!***
334 : :
335 : : end module m_gwls_wf
336 : : !!***
|
