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 : !!***
|