Line data Source code
1 : MODULE m_vacfun
2 : use m_juDFT
3 : #ifdef CPP_MPI
4 : use mpi
5 : #endif
6 : CONTAINS
7 288 : SUBROUTINE vacfun(&
8 : fmpi,vacuum,stars,input,nococonv,jspin1,jspin2,&
9 576 : cell,ivac,evac,bkpt, vxy,vz,kvac,nv2,&
10 576 : tuuv,tddv,tudv,tduv,uz,duz,udz,dudz,ddnv,wronk,&
11 864 : bkptq,v1,kvacq,nv2q,uzq,duzq,udzq,dudzq)
12 : !*********************************************************************
13 : ! determines the necessary values and derivatives on the vacuum
14 : ! boundary (ivac=1 upper vacuum; ivac=2, lower) for energy
15 : ! parameter evac. also sets up the 2d hamiltonian matrices
16 : ! necessary to update the full hamiltonian matrix.
17 : ! m. weinert
18 : !*********************************************************************
19 :
20 : USE m_constants
21 : USE m_intgr, ONLY : intgz0
22 : USE m_vacuz
23 : USE m_vacudz
24 : USE m_types
25 : IMPLICIT NONE
26 :
27 : TYPE(t_mpi),INTENT(IN) :: fmpi
28 : TYPE(t_input),INTENT(IN) :: input
29 : TYPE(t_vacuum),INTENT(IN) :: vacuum
30 : TYPE(t_nococonv),INTENT(IN) :: nococonv
31 : TYPE(t_stars),INTENT(IN) :: stars
32 : TYPE(t_cell),INTENT(IN) :: cell
33 : ! ..
34 : ! .. Scalar Arguments ..
35 : INTEGER, INTENT (IN) :: ivac,jspin1,jspin2
36 : REAL, INTENT (OUT) :: wronk
37 : ! ..
38 : ! .. Array Arguments ..
39 : INTEGER, INTENT (IN) :: nv2(:)!(input%jspins)
40 : INTEGER, INTENT (IN) :: kvac(:,:,:)!(2,lapw%dim_nv2d(),input%jspins)
41 : COMPLEX, INTENT (IN) :: vxy(:,:,:,:) !(vacuum%nmzxyd,stars%ng2-1,nvac,:)
42 : COMPLEX, INTENT (OUT):: tddv(:,:),tduv(:,:)!(lapw%dim_nv2d(),lapw%dim_nv2d())
43 : COMPLEX, INTENT (OUT):: tudv(:,:),tuuv(:,:)!(lapw%dim_nv2d(),lapw%dim_nv2d())
44 : COMPLEX, INTENT (IN) :: vz(:,:,:) !(vacuum%nmzd,2,4) ,
45 : REAL, INTENT (IN) :: evac(:,:)!(2,input%jspins)
46 : REAL, INTENT (IN) :: bkpt(3)
47 : REAL, INTENT (OUT):: udz(:,:),uz(:,:)!(lapw%dim_nv2d(),input%jspins)
48 : REAL, INTENT (OUT):: dudz(:,:),duz(:,:)!(lapw%dim_nv2d(),input%jspins)
49 : REAL, INTENT (OUT):: ddnv(:,:)!(lapw%dim_nv2d(),input%jspins)
50 : ! Optional DFPT stuff
51 : REAL, OPTIONAL, INTENT (IN) :: bkptq(3)
52 : REAL, OPTIONAL, INTENT (OUT):: udzq(:,:),uzq(:,:)!(lapw%dim_nv2d(),input%jspins)
53 : REAL, OPTIONAL, INTENT (OUT):: dudzq(:,:),duzq(:,:)!(lapw%dim_nv2d(),input%jspins)
54 : INTEGER, OPTIONAL, INTENT (IN) :: kvacq(:,:,:)!(2,lapw%dim_nv2d(),input%jspins)
55 : INTEGER, OPTIONAL, INTENT (IN) :: nv2q(:)!(input%jspins)
56 : COMPLEX, OPTIONAL, INTENT (IN) :: v1(:,:,:,:) !(vacuum%nmzxyd,stars%ng2-1,nvac,:)
57 : ! ..
58 : ! .. Local Scalars ..
59 : REAL ev,scale,xv,yv,vzero,fac
60 : COMPLEX phase
61 : INTEGER i,i1,i2,i3,ik,ind2,ind3,jk,np1,jspin,ipot,nbuf,ierr,loclen
62 : INTEGER mat_start,mat_end
63 : LOGICAL tail, l_dfpt
64 : ! ..
65 : ! .. Local Arrays ..
66 288 : REAL u(vacuum%nmzd,size(duz,1),input%jspins),ud(vacuum%nmzd,size(duz,1),input%jspins)
67 288 : REAL v(3),x(vacuum%nmzd), qssbti(2,2)
68 288 : COMPLEX, ALLOCATABLE :: tddv_loc(:,:), tduv_loc(:,:), tudv_loc(:,:), tuuv_loc(:,:)
69 288 : COMPLEX, ALLOCATABLE :: tv_gather_buf(:)
70 288 : REAL :: ddnvq(SIZE(ddnv,1),input%jspins)
71 : ! DFPT
72 288 : REAL uq(vacuum%nmzd,size(duz,1),input%jspins),udq(vacuum%nmzd,size(duz,1),input%jspins)
73 : ! ..
74 288 : l_dfpt = .FALSE.
75 288 : IF (PRESENT(bkptq)) l_dfpt = .TRUE.
76 288 : fac=MERGE(1.0,-1.0,jspin1>=jspin2)
77 576 : ipot=MERGE(jspin1,3,jspin1==jspin2)
78 :
79 5013952 : tuuv=0.0;tudv=0.0;tddv=0.0;tduv=0.0
80 128420 : udz=0.0;duz=0.0;ddnv=0.0;udz=0.;uz=0.
81 288 : tail = .true.
82 288 : IF (l_dfpt) THEN
83 0 : udzq=0.0;duzq=0.0;ddnvq=0.0;udzq=0.;uzq=0.
84 : END IF
85 288 : np1 = vacuum%nmzxy + 1
86 : !---> wronksian for the schrodinger equation given by an identity
87 288 : wronk = 2.0
88 : !---> setup the spin-spiral q-vector
89 864 : qssbti(1:2,1) = - nococonv%qss(1:2)/2
90 864 : qssbti(1:2,2) = + nococonv%qss(1:2)/2
91 : !---> generate basis functions for each 2-d k+g
92 576 : DO jspin = MIN(jspin1,jspin2),MAX(jspin1,jspin2)
93 18056 : DO ik = 1,nv2(jspin)
94 53304 : v(1:2) = bkpt(1:2) + kvac(:,ik,jspin) + qssbti(1:2,jspin)
95 17768 : v(3) = 0.0
96 302056 : ev = evac(ivac,jspin) - 0.5*dot_product(v,matmul(v,cell%bbmat))
97 17768 : vzero = vz(vacuum%nmzd,ivac,jspin)
98 : CALL vacuz(ev,REAL(vz(1:,ivac,jspin)),vzero,vacuum%nmz,vacuum%delz,&
99 4477536 : uz(ik,jspin),duz(ik,jspin),u(1,ik,jspin))
100 : CALL vacudz(ev,REAL(vz(1:,ivac,jspin)),vzero,vacuum%nmz,vacuum%delz,&
101 : udz(ik,jspin),dudz(ik,jspin),ddnv(ik,jspin),&
102 4459768 : ud(1,ik,jspin),duz(ik,jspin),u(1,ik,jspin))
103 : !---> make sure the solutions satisfy the wronksian
104 : scale = wronk/ (udz(ik,jspin)*duz(ik,jspin)-&
105 17768 : & dudz(ik,jspin)*uz(ik,jspin))
106 17768 : udz(ik,jspin) = scale*udz(ik,jspin)
107 17768 : dudz(ik,jspin) = scale*dudz(ik,jspin)
108 17768 : ddnv(ik,jspin) = scale*ddnv(ik,jspin)
109 4460056 : ud(:,ik,jspin) = scale*ud(:,ik,jspin)
110 : enddo
111 576 : if (l_dfpt) then
112 0 : DO ik = 1,nv2q(jspin)
113 0 : v(1:2) = bkptq(1:2) + kvacq(:,ik,jspin) + qssbti(1:2,jspin)
114 0 : v(3) = 0.0
115 0 : ev = evac(ivac,jspin) - 0.5*dot_product(v,matmul(v,cell%bbmat))
116 0 : vzero = vz(vacuum%nmzd,ivac,jspin)
117 : CALL vacuz(ev,REAL(vz(1:,ivac,jspin)),vzero,vacuum%nmz,vacuum%delz,&
118 0 : uzq(ik,jspin),duzq(ik,jspin),uq(1,ik,jspin))
119 : CALL vacudz(ev,REAL(vz(1:,ivac,jspin)),vzero,vacuum%nmz,vacuum%delz,&
120 : udzq(ik,jspin),dudzq(ik,jspin),ddnvq(ik,jspin),&
121 0 : udq(1,ik,jspin),duzq(ik,jspin),uq(1,ik,jspin))
122 : !---> make sure the solutions satisfy the wronksian
123 : scale = wronk/ (udzq(ik,jspin)*duzq(ik,jspin)-&
124 0 : & dudzq(ik,jspin)*uzq(ik,jspin)) !! TODO: Output wronsk always = 2.0?
125 0 : udzq(ik,jspin) = scale*udzq(ik,jspin)
126 0 : dudzq(ik,jspin) = scale*dudzq(ik,jspin)
127 0 : ddnvq(ik,jspin) = scale*ddnvq(ik,jspin)
128 0 : udq(:,ik,jspin) = scale*udq(:,ik,jspin)
129 : enddo
130 : end if
131 : ENDDO
132 288 : loclen = size(tddv,2)/fmpi%n_size + 1
133 288 : mat_start = fmpi%n_rank*loclen + 1
134 288 : mat_end = (fmpi%n_rank+1)*loclen
135 1152 : ALLOCATE(tddv_loc(size(tddv,1),mat_start:mat_end))
136 864 : ALLOCATE(tduv_loc(size(tddv,1),mat_start:mat_end))
137 864 : ALLOCATE(tudv_loc(size(tddv,1),mat_start:mat_end))
138 864 : ALLOCATE(tuuv_loc(size(tddv,1),mat_start:mat_end))
139 3247344 : tuuv_loc=0.0;tudv_loc=0.0;tddv_loc=0.0;tduv_loc=0.0
140 :
141 : !---> set up the tuuv, etc. matrices
142 12592 : DO jk = mat_start,mat_end
143 12520 : IF (jk>nv2(jspin2)) EXIT
144 12592 : IF (.NOT.l_dfpt) THEN
145 : !$OMP PARALLEL DO DEFAULT(none) &
146 : !$OMP& SHARED(tuuv_loc,tddv_loc,tudv_loc,tduv_loc,ddnv,vz,jk) &
147 : !$OMP& SHARED(stars,jspin1,jspin2,evac,nv2,kvac,vacuum,u,vxy,tail,fac,np1,ivac,ipot,ud) &
148 12304 : !$OMP& PRIVATE(i1,i2,i3,ind3,phase,ind2,x,xv,yv)
149 : DO ik = 1,nv2(jspin1)
150 :
151 : !---> determine the warping component of the potential
152 : i1 = fac*(kvac(1,ik,jspin1) - kvac(1,jk,jspin2))
153 : i2 = fac*(kvac(2,ik,jspin1) - kvac(2,jk,jspin2))
154 : i3 = 0
155 : ind3 = stars%ig(i1,i2,i3)
156 : IF (ind3.EQ.0) CYCLE
157 : phase = stars%rgphs(i1,i2,i3)
158 : ind2 = stars%ig2(ind3)
159 : IF (ind2.EQ.0) THEN
160 : WRITE (oUnit,FMT=8000) ik,jk
161 : 8000 FORMAT (' **** error in map2 for 2-d stars',2i5)
162 : CALL juDFT_error("error in map2 for 2-d stars",calledby ="vacfun")
163 : END IF
164 : !---> get the proper warping index (vxy starts with the 2nd star)
165 : ind2 = ind2 - 1
166 : IF (ind2.NE.0) THEN
167 : !---> only the warping part, 1st star (G=0) is done later
168 :
169 : !---> obtain the warping matrix elements
170 : !---> note that the tail correction (tail=.true.) is included for
171 : !---> the integrals, i.e. the integrand is from infinity inward
172 :
173 : !---> tuuv
174 : DO i = 1,vacuum%nmzxy
175 : x(np1-i) = u(i,ik,jspin1)*u(i,jk,jspin2)*REAL(vxy(i,ind2,ivac,ipot))
176 : enddo
177 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
178 : DO i = 1,vacuum%nmzxy
179 : x(np1-i) = u(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(vxy(i,ind2,ivac,ipot))
180 : enddo
181 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
182 : tuuv_loc(ik,jk) = phase*cmplx(xv,yv)
183 :
184 : !---> tddv
185 : DO i = 1,vacuum%nmzxy
186 : x(np1-i) = ud(i,ik,jspin1)*ud(i,jk,jspin2)*REAL(vxy(i,ind2,ivac,ipot))
187 : enddo
188 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
189 : DO i = 1,vacuum%nmzxy
190 : x(np1-i) =ud(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(vxy(i,ind2,ivac,ipot))
191 : enddo
192 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
193 : tddv_loc(ik,jk) = phase*cmplx(xv,yv)
194 :
195 : !---> tudv
196 : DO i = 1,vacuum%nmzxy
197 : x(np1-i) = u(i,ik,jspin1)*ud(i,jk,jspin2)*real(vxy(i,ind2,ivac,ipot))
198 : enddo
199 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
200 : DO i = 1,vacuum%nmzxy
201 : x(np1-i) = u(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(vxy(i,ind2,ivac,ipot))
202 : enddo
203 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
204 : tudv_loc(ik,jk) = phase*cmplx(xv,yv)
205 :
206 : !---> tduv
207 : DO i = 1,vacuum%nmzxy
208 : x(np1-i) = ud(i,ik,jspin1)*u(i,jk,jspin2)*real(vxy(i,ind2,ivac,ipot))
209 : enddo
210 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
211 : DO i = 1,vacuum%nmzxy
212 : x(np1-i) = ud(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(vxy(i,ind2,ivac,ipot))
213 : enddo
214 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
215 : tduv_loc(ik,jk) = phase*cmplx(xv,yv)
216 :
217 : ELSE
218 : !---> diagonal (film muffin-tin) terms
219 : IF (jspin1==jspin2) THEN
220 : tuuv_loc(ik,ik) = cmplx(evac(ivac,jspin1),0.0)
221 : tddv_loc(ik,ik) = cmplx(evac(ivac,jspin1)*ddnv(ik,jspin1),0.0)
222 : tudv_loc(ik,ik) = cmplx(0.5,0.0)
223 : tduv_loc(ik,ik) = cmplx(0.5,0.0)
224 : ELSE
225 :
226 : !---> tuuv
227 : DO i = 1,vacuum%nmz
228 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*u(i,jk,jspin2)*vz(i,ivac,3)
229 : ENDDO
230 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
231 : DO i = 1,vacuum%nmz
232 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(vz(i,ivac,3))
233 : ENDDO
234 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
235 : tuuv_loc(ik,jk) = cmplx(xv,yv)
236 :
237 : !---> tddv
238 : DO i = 1,vacuum%nmz
239 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*ud(i,jk,jspin2)*vz(i,ivac,3)
240 : ENDDO
241 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
242 : DO i = 1,vacuum%nmz
243 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(vz(i,ivac,3))
244 : ENDDO
245 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
246 : tddv_loc(ik,jk) = cmplx(xv,yv)
247 :
248 : !---> tudv
249 : DO i = 1,vacuum%nmz
250 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*ud(i,jk,jspin2)*vz(i,ivac,3)
251 : ENDDO
252 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
253 : DO i = 1,vacuum%nmz
254 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(vz(i,ivac,3))
255 : ENDDO
256 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
257 : tudv_loc(ik,jk) = cmplx(xv,yv)
258 :
259 : !---> tduv
260 : DO i = 1,vacuum%nmz
261 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*u(i,jk,jspin2)*vz(i,ivac,3)
262 : ENDDO
263 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
264 : DO i = 1,vacuum%nmz
265 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(vz(i,ivac,3))
266 : ENDDO
267 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
268 : tduv_loc(ik,jk) = cmplx(xv,yv)
269 : ENDIF
270 :
271 : ENDIF
272 : ENDDO
273 : !$OMP END PARALLEL DO
274 : ELSE
275 : !$OMP PARALLEL DO DEFAULT(none) &
276 : !$OMP& SHARED(tuuv_loc,tddv_loc,tudv_loc,tduv_loc,ddnv,vz,v1,jk,bkpt,bkptq) &
277 : !$OMP& SHARED(stars,jspin1,jspin2,evac,nv2,nv2q,kvac,kvacq,vacuum,u,uq,tail,fac,np1,ivac,ipot,ud,udq,l_dfpt) &
278 0 : !$OMP& PRIVATE(i1,i2,i3,ind3,phase,ind2,x,xv,yv)
279 : DO ik = 1,nv2q(jspin1)
280 : !---> determine the warping component of the potential
281 : i1 = fac*(kvacq(1,ik,jspin1) - kvac(1,jk,jspin2))
282 : i2 = fac*(kvacq(2,ik,jspin1) - kvac(2,jk,jspin2))
283 : i3 = 0
284 : ind3 = stars%ig(i1,i2,i3)
285 : IF (ind3.EQ.0) CYCLE
286 : phase = stars%rgphs(i1,i2,i3)
287 : ind2 = stars%ig2(ind3)
288 : IF (ind2.EQ.0) THEN
289 : WRITE (oUnit,FMT=8001) ik,jk
290 : 8001 FORMAT (' **** error in map2 for 2-d stars',2i5)
291 : CALL juDFT_error("error in map2 for 2-d stars",calledby ="vacfun")
292 : END IF
293 : !---> get the proper warping index (v1xy starts with the 2nd star)
294 : ind2 = ind2 - 1
295 : IF ((ind2.NE.0).OR.(norm2(bkptq-bkpt)>1e-8)) THEN
296 : !---> tuuv
297 : DO i = 1,vacuum%nmzxy
298 : x(np1-i) = uq(i,ik,jspin1)*u(i,jk,jspin2)*REAL(v1(i,ind2+1,ivac,ipot))
299 : enddo
300 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
301 : DO i = 1,vacuum%nmzxy
302 : x(np1-i) = uq(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(v1(i,ind2+1,ivac,ipot))
303 : enddo
304 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
305 : tuuv_loc(ik,jk) = phase*cmplx(xv,yv)
306 :
307 : !---> tddv
308 : DO i = 1,vacuum%nmzxy
309 : x(np1-i) = udq(i,ik,jspin1)*ud(i,jk,jspin2)*REAL(v1(i,ind2+1,ivac,ipot))
310 : enddo
311 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
312 : DO i = 1,vacuum%nmzxy
313 : x(np1-i) =udq(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(v1(i,ind2+1,ivac,ipot))
314 : enddo
315 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
316 : tddv_loc(ik,jk) = phase*cmplx(xv,yv)
317 :
318 : !---> tudv
319 : DO i = 1,vacuum%nmzxy
320 : x(np1-i) = uq(i,ik,jspin1)*ud(i,jk,jspin2)*real(v1(i,ind2+1,ivac,ipot))
321 : enddo
322 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
323 : DO i = 1,vacuum%nmzxy
324 : x(np1-i) = uq(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(v1(i,ind2+1,ivac,ipot))
325 : enddo
326 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
327 : tudv_loc(ik,jk) = phase*cmplx(xv,yv)
328 :
329 : !---> tduv
330 : DO i = 1,vacuum%nmzxy
331 : x(np1-i) = udq(i,ik,jspin1)*u(i,jk,jspin2)*real(v1(i,ind2+1,ivac,ipot))
332 : enddo
333 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,xv,tail)
334 : DO i = 1,vacuum%nmzxy
335 : x(np1-i) = udq(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(v1(i,ind2+1,ivac,ipot))
336 : enddo
337 : CALL intgz0(x,vacuum%delz,vacuum%nmzxy,yv,tail)
338 : tduv_loc(ik,jk) = phase*cmplx(xv,yv)
339 :
340 : ELSE
341 : !---> diagonal (film muffin-tin) terms
342 : !---> tuuv
343 : DO i = 1,vacuum%nmz
344 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*u(i,jk,jspin2)*v1(i,1,ivac,ipot)
345 : ENDDO
346 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
347 : DO i = 1,vacuum%nmz
348 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(v1(i,1,ivac,ipot))
349 : ENDDO
350 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
351 : tuuv_loc(ik,jk) = cmplx(xv,yv)
352 :
353 : !---> tddv
354 : DO i = 1,vacuum%nmz
355 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*ud(i,jk,jspin2)*v1(i,1,ivac,ipot)
356 : ENDDO
357 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
358 : DO i = 1,vacuum%nmz
359 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(v1(i,1,ivac,ipot))
360 : ENDDO
361 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
362 : tddv_loc(ik,jk) = cmplx(xv,yv)
363 :
364 : !---> tudv
365 : DO i = 1,vacuum%nmz
366 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*ud(i,jk,jspin2)*v1(i,1,ivac,ipot)
367 : ENDDO
368 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
369 : DO i = 1,vacuum%nmz
370 : x(vacuum%nmz+1-i) = u(i,ik,jspin1)*ud(i,jk,jspin2)*fac*AIMAG(v1(i,1,ivac,ipot))
371 : ENDDO
372 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
373 : tudv_loc(ik,jk) = cmplx(xv,yv)
374 :
375 : !---> tduv
376 : DO i = 1,vacuum%nmz
377 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*u(i,jk,jspin2)*v1(i,1,ivac,ipot)
378 : ENDDO
379 : CALL intgz0(x,vacuum%delz,vacuum%nmz,xv,tail)
380 : DO i = 1,vacuum%nmz
381 : x(vacuum%nmz+1-i) = ud(i,ik,jspin1)*u(i,jk,jspin2)*fac*AIMAG(v1(i,1,ivac,ipot))
382 : ENDDO
383 : CALL intgz0(x,vacuum%delz,vacuum%nmz,yv,tail)
384 : tduv_loc(ik,jk) = cmplx(xv,yv)
385 : ENDIF
386 : ENDDO
387 : !$OMP END PARALLEL DO
388 : END IF
389 : ENDDO
390 :
391 288 : IF ( fmpi%n_size == 1 ) THEN
392 344560 : tuuv = tuuv_loc(:,:size(tuuv,2))
393 344560 : tduv = tduv_loc(:,:size(tduv,2))
394 344560 : tudv = tudv_loc(:,:size(tudv,2))
395 344560 : tddv = tddv_loc(:,:size(tddv,2))
396 : ELSE
397 432 : nbuf = size(tuuv_loc)
398 : #ifdef CPP_MPI
399 432 : ALLOCATE(tv_gather_buf(nbuf*fmpi%n_size))
400 144 : CALL MPI_ALLGATHER(tuuv_loc,nbuf,MPI_DOUBLE_COMPLEX,tv_gather_buf,nbuf,MPI_DOUBLE_COMPLEX,fmpi%sub_comm,ierr)
401 432 : CALL zcopy(size(tuuv),tv_gather_buf,1,tuuv,1)
402 144 : CALL MPI_ALLGATHER(tduv_loc,nbuf,MPI_DOUBLE_COMPLEX,tv_gather_buf,nbuf,MPI_DOUBLE_COMPLEX,fmpi%sub_comm,ierr)
403 432 : CALL zcopy(size(tduv),tv_gather_buf,1,tduv,1)
404 144 : CALL MPI_ALLGATHER(tudv_loc,nbuf,MPI_DOUBLE_COMPLEX,tv_gather_buf,nbuf,MPI_DOUBLE_COMPLEX,fmpi%sub_comm,ierr)
405 432 : CALL zcopy(size(tudv),tv_gather_buf,1,tudv,1)
406 144 : CALL MPI_ALLGATHER(tddv_loc,nbuf,MPI_DOUBLE_COMPLEX,tv_gather_buf,nbuf,MPI_DOUBLE_COMPLEX,fmpi%sub_comm,ierr)
407 432 : CALL zcopy(size(tddv),tv_gather_buf,1,tddv,1)
408 144 : DEALLOCATE (tv_gather_buf)
409 : #endif
410 : ENDIF
411 :
412 288 : DEALLOCATE(tddv_loc, tduv_loc, tudv_loc, tuuv_loc)
413 :
414 288 : END SUBROUTINE vacfun
415 : END MODULE m_vacfun
|
