Line data Source code
1 : ! -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- --
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_hsmt_sph
8 : USE m_juDFT
9 : IMPLICIT NONE
10 :
11 : INTERFACE hsmt_sph
12 : #ifdef _OPENACC
13 : MODULE PROCEDURE hsmt_sph_acc
14 : #else
15 : MODULE PROCEDURE hsmt_sph_cpu
16 : #endif
17 : END INTERFACE
18 :
19 : CONTAINS
20 :
21 0 : SUBROUTINE hsmt_sph_acc(n,atoms,fmpi,isp,input,nococonv,igSpinPr,igSpin,chi,lapw,el,e_shift,usdus,fjgj,smat,hmat,set0,l_fullj,lapwq,fjgjq)
22 : USE m_constants, ONLY : fpi_const, tpi_const
23 : USE m_types
24 : USE m_hsmt_fjgj
25 :
26 :
27 : TYPE(t_input), INTENT(IN) :: input
28 : TYPE(t_mpi), INTENT(IN) :: fmpi
29 : TYPE(t_nococonv), INTENT(IN) :: nococonv
30 : TYPE(t_atoms), INTENT(IN) :: atoms
31 : TYPE(t_lapw), INTENT(IN) :: lapw
32 : TYPE(t_usdus), INTENT(IN) :: usdus
33 : TYPE(t_fjgj), INTENT(IN) :: fjgj
34 : CLASS(t_mat), INTENT(INOUT) :: smat, hmat
35 : LOGICAL, INTENT(IN) :: l_fullj, set0 !if true, initialize the smat matrix with zeros
36 :
37 : TYPE(t_lapw), OPTIONAL, INTENT(IN) :: lapwq
38 : TYPE(t_fjgj), OPTIONAL, INTENT(IN) :: fjgjq
39 :
40 : ! Scalar Arguments
41 : INTEGER, INTENT(IN) :: n, isp, igSpinPr, igSpin
42 : COMPLEX, INTENT(IN) :: chi
43 :
44 : ! Array Arguments
45 : REAL, INTENT(IN) :: el(0:atoms%lmaxd,atoms%ntype,input%jspins)
46 : REAL, INTENT(IN) :: e_shift!(atoms%ntype,input%jspins)
47 :
48 : ! Local Scalars
49 : REAL :: tnn(3), elall, fjkiln, gjkiln, ddnln, ski(3)
50 : REAL :: apw_lo1, apw_lo2, w1
51 :
52 : INTEGER :: ikG0, ikG, ikGPr, l, nn, l3, jv, kj_off, kj_vec
53 :
54 : ! Local Arrays
55 0 : REAL :: fleg1(0:atoms%lmaxd), fleg2(0:atoms%lmaxd), fl2p1(0:atoms%lmaxd)
56 : REAL :: qssAdd(3), qssAddPr(3)
57 : REAL :: plegend(0:2)
58 : REAL :: xlegend
59 : REAL :: VecHelpS, VecHelpH
60 : REAL :: cph_re, cph_im
61 : REAL :: dot, fct, fct2
62 :
63 : COMPLEX :: cfac
64 :
65 : LOGICAL :: l_samelapw
66 :
67 0 : TYPE(t_lapw) :: lapwPr
68 :
69 0 : CALL timestart("spherical setup")
70 0 : l_samelapw = .FALSE.
71 0 : IF (.NOT.PRESENT(lapwq)) l_samelapw = .TRUE.
72 : IF (.NOT.l_samelapw) THEN
73 0 : lapwPr = lapwq
74 : ELSE
75 0 : lapwPr = lapw
76 : END IF
77 : !call nvtxStartRange("hsmt_sph",1)
78 0 : DO l = 0,atoms%lmaxd
79 0 : fleg1(l) = REAL(l+l+1)/REAL(l+1)
80 0 : fleg2(l) = REAL(l)/REAL(l+1)
81 0 : fl2p1(l) = REAL(l+l+1)/fpi_const
82 : END DO ! l
83 :
84 0 : qssAdd = MERGE(-nococonv%qss/2, +nococonv%qss/2, igSpin.EQ.1)
85 0 : qssAddPr = MERGE(-nococonv%qss/2, +nococonv%qss/2, igSpinPr.EQ.1)
86 : !$acc data &
87 : !$acc& copyin(igSpin,igSpinPr,n,fleg1,fleg2,isp,fl2p1,el,e_shift,chi,qssAdd,qssAddPr,l_fullj)&
88 : !$acc& copyin(lapw,lapwPr,atoms,fmpi,input,usdus)&
89 : !$acc& copyin(lapw%nv,lapw%gvec,lapw%gk,lapwPr%nv,lapwPr%gvec,lapwPr%gk,lapw%bkpt,lapwPr%bkpt)&
90 : !$acc& copyin(atoms%lmax,atoms%rmt,atoms%lnonsph,atoms%firstAtom,atoms%neq,atoms%taual)&
91 : !$acc& copyin(fmpi%n_size,fmpi%n_rank)&
92 : !$acc& copyin(input%l_useapw)&
93 : !$acc& copyin(usdus%dus,usdus%uds,usdus%us,usdus%ddn,usdus%duds)&
94 : !$acc& present(fjgj)&
95 : !$acc& present(hmat,smat,hmat%data_c,hmat%data_r,smat%data_r,smat%data_c)
96 :
97 : !$acc parallel default(none)
98 : !$acc loop gang
99 0 : DO ikG = fmpi%n_rank+1, lapw%nv(igSpin), fmpi%n_size
100 : !$acc loop vector independent&
101 : !$acc & PRIVATE(ikGPr,ikG0,ski,plegend,tnn,vechelps,vechelph,xlegend,fjkiln,gjkiln,ddnln,elall,l3,l,fct,fct2,cph_re,cph_im,cfac,dot)
102 0 : DO ikGPr = 1, MERGE(lapwPr%nv(igSpinPr),MIN(ikG,lapwPr%nv(igSpinPr)),l_fullj)
103 0 : ikG0 = (ikG-1)/fmpi%n_size + 1
104 0 : ski = lapw%gvec(:,ikG,igSpin) + qssAdd(:) + lapw%bkpt + lapw%qphon
105 :
106 : ! Update overlap and l-diagonal hamiltonian matrix
107 0 : VecHelpS = 0.0
108 0 : VecHelpH = 0.0
109 :
110 : ! x for legendre polynomials
111 0 : xlegend = dot_product(lapwPr%gk(1:3,ikGPr,igSpinPr),lapw%gk(1:3,ikG,igSpin))
112 :
113 : !$acc loop seq
114 0 : DO l = 0,atoms%lmax(n)
115 0 : fjkiln = fjgj%fj(ikG,l,isp,igSpin)
116 0 : gjkiln = fjgj%gj(ikG,l,isp,igSpin)
117 0 : IF (input%l_useapw) THEN
118 0 : w1 = 0.5 * ( usdus%uds(l,n,isp)*usdus%dus(l,n,isp) + usdus%us(l,n,isp)*usdus%duds(l,n,isp) )
119 0 : apw_lo1 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( gjkiln * w1 + fjkiln * usdus%us(l,n,isp) * usdus%dus(l,n,isp) )
120 0 : apw_lo2 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( fjkiln * w1 + gjkiln * usdus%uds(l,n,isp) * usdus%duds(l,n,isp) )
121 : END IF ! useapw
122 0 : ddnln = usdus%ddn(l,n,isp)
123 0 : elall = el(l,n,isp)
124 :
125 0 : IF (l<=atoms%lnonsph(n).AND..NOT.l_fullj) elall=elall-e_shift!(isp)
126 :
127 : ! Legendre polynomials
128 0 : l3 = modulo(l, 3)
129 :
130 0 : IF (l == 0) THEN
131 0 : plegend(0) = 1.0
132 0 : ELSE IF (l == 1) THEN
133 0 : plegend(1) = xlegend
134 : ELSE
135 : plegend(l3) = fleg1(l-1)*xlegend*plegend(modulo(l-1,3)) &
136 0 : & - fleg2(l-1)*plegend(modulo(l-2,3))
137 : END IF ! l
138 :
139 : fct = plegend(l3) * fl2p1(l) * ( fjkiln*fjgj%fj(ikGPr,l,isp,igSpinPr) &
140 0 : & + gjkiln*fjgj%gj(ikGPr,l,isp,igSpinPr)*ddnln )
141 : !IF (.NOT.l_fullj) THEN
142 0 : IF (.TRUE.) THEN
143 : fct2 = plegend(l3)*fl2p1(l) * 0.5 * ( gjkiln*fjgj%fj(ikGPr,l,isp,igSpinPr) &
144 0 : & + fjkiln*fjgj%gj(ikGPr,l,isp,igSpinPr) )
145 : ELSE
146 : fct2 = plegend(l3)*fl2p1(l) * gjkiln*fjgj%fj(ikGPr,l,isp,igSpinPr)
147 : END IF
148 :
149 0 : VecHelpS = VecHelpS + fct
150 0 : VecHelpH = VecHelpH + fct*elall + fct2
151 :
152 0 : IF (input%l_useapw) THEN
153 : VecHelpH = VecHelpH + plegend(l3) * ( apw_lo1*fjgj%fj(ikGPr,l,isp,igSpinPr) &
154 0 : & + apw_lo2*fjgj%gj(l,ikGPr,isp,igSpinPr) )
155 : END IF ! useapw
156 : END DO ! l
157 : !$end acc
158 : ! Set up phase factors
159 0 : cph_re = 0.0
160 0 : cph_im = 0.0
161 0 : DO nn = atoms%firstAtom(n), atoms%firstAtom(n) + atoms%neq(n) - 1
162 0 : tnn(1:3) = tpi_const*atoms%taual(1:3,nn)
163 :
164 0 : dot = DOT_PRODUCT(ski(1:3) - lapwPr%gvec(1:3,ikGPr,igSpinPr) - qssAddPr(1:3) - lapwPr%bkpt - lapwPr%qphon, tnn(1:3))
165 :
166 0 : cph_re = cph_re + COS(dot)
167 0 : cph_im = cph_im + SIN(dot)
168 : ! IF (igSpinPr.NE.igSpin) cph_im=-cph_im
169 : END DO ! nn
170 :
171 0 : cfac = CMPLX(cph_re,cph_im)
172 : ! ! Prefactor: i * (k + G + qssAdd - k' - G' - qssAdd')
173 : ! IF (l_fullj) THEN
174 : ! pref = ImagUnit * MATMUL(ski(1:3) - lapwPr%gvec(1:3,ikGPr,igSpinPr) - qssAddPr(1:3) - lapwPr%bkpt, bmat)
175 : ! cfac = pref(idir) * cfac
176 : ! END IF
177 :
178 0 : IF (smat%l_real) THEN
179 0 : IF (set0) THEN
180 0 : smat%data_r(ikGPr,ikG0) = cph_re * VecHelpS
181 : ELSE
182 : smat%data_r(ikGPr,ikG0) = &
183 0 : smat%data_r(ikGPr,ikG0) + cph_re * VecHelpS
184 : END IF
185 : hmat%data_r(ikGPr,ikG0) = &
186 0 : hmat%data_r(ikGPr,ikG0) + cph_re * VecHelpH
187 : ELSE ! real
188 0 : IF (set0) THEN
189 0 : smat%data_c(ikGPr,ikG0) = chi*cfac * VecHelpS
190 : ELSE
191 : smat%data_c(ikGPr,ikG0) = &
192 0 : smat%data_c(ikGPr,ikG0) + chi*cfac * VecHelpS
193 : END IF
194 : hmat%data_c(ikGPr,ikG0) = &
195 0 : hmat%data_c(ikGPr,ikG0) + chi*cfac * VecHelpH
196 : END IF ! real
197 : END DO ! kj_off
198 : !$acc end loop
199 : END DO ! ikG
200 :
201 : !$acc end loop
202 : !$acc end parallel
203 : !$acc end data
204 : !$acc wait
205 0 : CALL timestop("spherical setup")
206 : !call nvtxEndRange()
207 0 : RETURN
208 0 : END SUBROUTINE hsmt_sph_acc
209 :
210 17504 : SUBROUTINE hsmt_sph_cpu(n,atoms,fmpi,isp,input,nococonv,igSpinPr,igSpin,chi,lapw,el,e_shift,usdus,fjgj,smat,hmat,set0,l_fullj,lapwq, fjgjq)
211 : USE m_constants, ONLY : fpi_const, tpi_const
212 : USE m_types
213 : USE m_hsmt_fjgj
214 :
215 :
216 : TYPE(t_input), INTENT(IN) :: input
217 : TYPE(t_mpi), INTENT(IN) :: fmpi
218 : TYPE(t_nococonv), INTENT(IN) :: nococonv
219 : TYPE(t_atoms), INTENT(IN) :: atoms
220 : TYPE(t_lapw), INTENT(IN) :: lapw
221 : TYPE(t_usdus), INTENT(IN) :: usdus
222 : TYPE(t_fjgj), INTENT(IN) :: fjgj
223 : CLASS(t_mat), INTENT(INOUT) :: smat, hmat
224 : LOGICAL, INTENT(IN) :: l_fullj, set0 !if true, initialize the smat matrix with zeros
225 :
226 : TYPE(t_lapw), OPTIONAL, INTENT(IN) :: lapwq
227 : TYPE(t_fjgj), OPTIONAL, INTENT(IN) :: fjgjq
228 :
229 : ! Scalar Arguments
230 : INTEGER, INTENT(IN) :: n, isp, igSpinPr, igSpin
231 : COMPLEX, INTENT(IN) :: chi
232 :
233 : ! Array Arguments
234 : REAL, INTENT(IN) :: el(0:atoms%lmaxd,atoms%ntype,input%jspins)
235 : REAL, INTENT(IN) :: e_shift!(atoms%ntype,input%jspins)
236 :
237 : ! Local Scalars
238 : REAL :: tnn(3), elall, fjkiln, gjkiln, ddnln, ski(3)
239 : REAL :: apw_lo1, apw_lo2, w1
240 :
241 : INTEGER :: ikG0, ikG, ikGPr, l, nn, kj_end, l3, jv, kj_off, kj_vec
242 :
243 : LOGICAL :: l_samelapw
244 :
245 17504 : TYPE(t_lapw) :: lapwPr
246 17504 : TYPE(t_fjgj) :: fjgjPr
247 :
248 : ! Local Arrays
249 17504 : REAL :: fleg1(0:atoms%lmaxd),fleg2(0:atoms%lmaxd),fl2p1(0:atoms%lmaxd)
250 : REAL :: qssAdd(3),qssAddPr(3)
251 :
252 17504 : REAL, ALLOCATABLE :: plegend(:,:)
253 17504 : REAL, ALLOCATABLE :: xlegend(:)
254 17504 : REAL, ALLOCATABLE :: VecHelpS(:),VecHelpH(:)
255 17504 : REAL, ALLOCATABLE :: cph_re(:), cph_im(:)
256 17504 : REAL, ALLOCATABLE :: dot(:), fct(:), fct2(:)
257 :
258 17504 : COMPLEX, ALLOCATABLE :: cfac(:)
259 :
260 : INTEGER, PARAMETER :: NVEC = 128
261 :
262 : INTEGER :: NVEC_rem !remainder
263 :
264 17504 : CALL timestart("spherical setup")
265 17504 : l_samelapw = .FALSE.
266 17504 : IF (.NOT.PRESENT(lapwq)) l_samelapw = .TRUE.
267 : IF (.NOT.l_samelapw) THEN
268 0 : lapwPr = lapwq
269 0 : fjgjPr = fjgjq
270 : ELSE
271 17504 : lapwPr = lapw
272 17504 : fjgjPr = fjgj
273 : END IF
274 : !call nvtxStartRange("hsmt_sph",1)
275 193876 : DO l = 0,atoms%lmaxd
276 176372 : fleg1(l) = REAL(l+l+1)/REAL(l+1)
277 176372 : fleg2(l) = REAL(l)/REAL(l+1)
278 193876 : fl2p1(l) = REAL(l+l+1)/fpi_const
279 : END DO ! l
280 :
281 : !$OMP PARALLEL DEFAULT(NONE)&
282 : !$OMP SHARED(lapw,lapwPr,atoms,nococonv,fmpi,input,usdus,smat,hmat)&
283 : !$OMP SHARED(igSpin,igSpinPr,n,fleg1,fleg2,fjgj,fjgjPr,isp,fl2p1,el,e_shift,chi,set0,l_fullj)&
284 : !$OMP PRIVATE(ikG0,ikG,ski,ikGPr,kj_off,kj_vec,plegend,xlegend,l,l3,kj_end,qssAdd,qssAddPr,fct2)&
285 : !$OMP PRIVATE(cph_re,cph_im,cfac,dot,nn,tnn,fjkiln,gjkiln)&
286 : !$OMP PRIVATE(w1,apw_lo1,apw_lo2,ddnln,elall,fct)&
287 17504 : !$OMP PRIVATE(VecHelpS,VecHelpH,NVEC_rem)
288 : ALLOCATE(cph_re(NVEC),cph_im(NVEC),cfac(NVEC))
289 : ALLOCATE(dot(NVEC),fct(NVEC),fct2(NVEC))
290 : ALLOCATE(plegend(NVEC,0:2))
291 : ALLOCATE(xlegend(NVEC))
292 : ALLOCATE(VecHelpS(NVEC),VecHelpH(NVEC))
293 : qssAdd = MERGE(-nococonv%qss/2,+nococonv%qss/2,igSpin .EQ.1)
294 : qssAddPr = MERGE(-nococonv%qss/2,+nococonv%qss/2,igSpinPr.EQ.1)
295 : !$OMP DO SCHEDULE(DYNAMIC,1)
296 : DO ikG = fmpi%n_rank+1, lapw%nv(igSpin), fmpi%n_size
297 : kj_end = MERGE(lapwPr%nv(igSpinPr),min(ikG,lapwPr%nv(igSpinPr)),l_fullj)
298 : ikG0 = (ikG-1)/fmpi%n_size + 1
299 : ski = lapw%gvec(:,ikG,igSpin) + qssAdd(:) + lapw%bkpt + lapw%qphon
300 : DO kj_off = 1, lapwPr%nv(igSpinPr), NVEC
301 : NVEC_rem = NVEC
302 : kj_vec = kj_off - 1 + NVEC
303 : IF (kj_vec > kj_end) THEN
304 : kj_vec = kj_end
305 : NVEC_rem = kj_end - kj_off + 1
306 : END IF
307 : IF (NVEC_rem<0 ) exit
308 : ! Update overlap and l-diagonal hamiltonian matrix
309 : VecHelpS = 0.0
310 : VecHelpH = 0.0
311 :
312 : ! x for legendre polynomials
313 : DO jv = 0, NVEC_rem-1
314 : ikGPr = jv + kj_off
315 : xlegend(jv+1) =dot_product(lapwPr%gk(1:3,ikGPr,igSpinPr),lapw%gk(1:3,ikG,igSpin))
316 : END DO ! ikGPr
317 :
318 : DO l = 0,atoms%lmax(n)
319 : fjkiln = fjgj%fj(ikG,l,isp,igSpin)
320 : gjkiln = fjgj%gj(ikG,l,isp,igSpin)
321 :
322 : IF (input%l_useapw) THEN
323 : w1 = 0.5 * ( usdus%uds(l,n,isp)*usdus%dus(l,n,isp) + usdus%us(l,n,isp)*usdus%duds(l,n,isp) )
324 : apw_lo1 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( gjkiln * w1 + fjkiln * usdus%us(l,n,isp) * usdus%dus(l,n,isp) )
325 : apw_lo2 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( fjkiln * w1 + gjkiln * usdus%uds(l,n,isp) * usdus%duds(l,n,isp) )
326 : END IF ! useapw
327 :
328 : !IF (l_fullj) THEN
329 : ! !apw_lo1 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( usdus%us(l,n,isp) * usdus%duds(l,n,isp) * gjkiln &
330 : ! ! & + usdus%us(l,n,isp) * usdus%dus(l,n,isp) * fjkiln )
331 : ! !apw_lo2 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( usdus%uds(l,n,isp) * usdus%dus(l,n,isp) * fjkiln &
332 : ! ! & + usdus%uds(l,n,isp) * usdus%duds(l,n,isp) * gjkiln)
333 : ! w1 = 0.5 * ( usdus%uds(l,n,isp)*usdus%dus(l,n,isp) + usdus%us(l,n,isp)*usdus%duds(l,n,isp) )
334 : ! apw_lo1 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( gjkiln * w1 + fjkiln * usdus%us(l,n,isp) * usdus%dus(l,n,isp) )
335 : ! apw_lo2 = fl2p1(l) * 0.5 * atoms%rmt(n)**2 * ( fjkiln * w1 + gjkiln * usdus%uds(l,n,isp) * usdus%duds(l,n,isp) )
336 : !END IF
337 :
338 : ddnln = usdus%ddn(l,n,isp)
339 : elall = el(l,n,isp)
340 :
341 : IF (l<=atoms%lnonsph(n).AND..NOT.l_fullj) elall=elall-e_shift!(isp)
342 : !IF (l<=atoms%lnonsph(n)) elall=elall-e_shift!(isp)
343 :
344 : ! Legendre polynomials
345 : l3 = modulo(l, 3)
346 : IF (l == 0) THEN
347 : plegend(:NVEC_REM,0) = 1.0
348 : ELSE IF (l == 1) THEN
349 : plegend(:NVEC_REM,1) = xlegend(:NVEC_REM)
350 : ELSE
351 : plegend(:NVEC_REM,l3) = fleg1(l-1)*xlegend(:NVEC_REM)*plegend(:NVEC_REM,modulo(l-1,3)) &
352 : & - fleg2(l-1)*plegend(:NVEC_REM,modulo(l-2,3))
353 : END IF ! l
354 :
355 : fct(:NVEC_REM) = plegend(:NVEC_REM,l3) * fl2p1(l) * ( fjkiln*fjgjPr%fj(kj_off:kj_vec,l,isp,igSpinPr) &
356 : & + gjkiln*fjgjPr%gj(kj_off:kj_vec,l,isp,igSpinPr)*ddnln )
357 :
358 : !IF (.NOT.l_fullj) THEN
359 : IF (.TRUE.) THEN
360 : fct2(:NVEC_REM) = plegend(:NVEC_REM,l3) * fl2p1(l) * 0.5 * ( gjkiln*fjgjPr%fj(kj_off:kj_vec,l,isp,igSpinPr) &
361 : & + fjkiln*fjgjPr%gj(kj_off:kj_vec,l,isp,igSpinPr) )
362 : ELSE
363 : fct2(:NVEC_REM) = plegend(:NVEC_REM,l3) * fl2p1(l) * gjkiln*fjgjPr%fj(kj_off:kj_vec,l,isp,igSpinPr)
364 : END IF
365 :
366 : VecHelpS(:NVEC_REM) = VecHelpS(:NVEC_REM) + fct(:NVEC_REM)
367 : VecHelpH(:NVEC_REM) = VecHelpH(:NVEC_REM) + fct(:NVEC_REM)*elall + fct2(:NVEC_REM)
368 :
369 : !IF (input%l_useapw.OR.(l_fullj.AND.l==0)) THEN
370 : !IF (input%l_useapw.OR.l_fullj) THEN ! correction
371 : IF (input%l_useapw) THEN
372 : VecHelpH(:NVEC_REM) = VecHelpH(:NVEC_REM) + plegend(:NVEC_REM,l3) * ( apw_lo1*fjgjPr%fj(kj_off:kj_vec,l,isp,igSpinPr) + &
373 : & apw_lo2*fjgjPr%gj(kj_off:kj_vec,l,isp,igSpinPr) )
374 : END IF ! useapw
375 : END DO ! l
376 : ! Set up phase factors
377 : cph_re = 0.0
378 : cph_im = 0.0
379 : DO nn = atoms%firstAtom(n), atoms%firstAtom(n) + atoms%neq(n) - 1
380 : tnn(1:3) = tpi_const*atoms%taual(1:3,nn)
381 : DO jv = 0, NVEC_rem-1
382 : ikGPr = jv + kj_off
383 : dot(jv+1) = DOT_PRODUCT(ski(1:3) - lapwPr%gvec(1:3,ikGPr,igSpinPr) - qssAddPr(1:3) - lapwPr%bkpt - lapwPr%qphon, tnn(1:3))
384 : END DO ! ikGPr
385 : cph_re(:NVEC_REM) = cph_re(:NVEC_REM) + COS(dot(:NVEC_REM))
386 : cph_im(:NVEC_REM) = cph_im(:NVEC_REM) + SIN(dot(:NVEC_REM))
387 : cfac(:NVEC_REM) = CMPLX(cph_re(:NVEC_REM),cph_im(:NVEC_REM))
388 : ! IF (igSpinPr.NE.igSpin) cph_im=-cph_im
389 : END DO ! nn
390 :
391 : IF (smat%l_real) THEN
392 : IF (set0) THEN
393 : smat%data_r(kj_off:kj_vec,ikG0) = cph_re(:NVEC_REM) * VecHelpS(:NVEC_REM)
394 : ELSE
395 : smat%data_r(kj_off:kj_vec,ikG0) = &
396 : smat%data_r(kj_off:kj_vec,ikG0) + cph_re(:NVEC_REM) * VecHelpS(:NVEC_REM)
397 : END IF
398 : hmat%data_r(kj_off:kj_vec,ikG0) = &
399 : hmat%data_r(kj_off:kj_vec,ikG0) + cph_re(:NVEC_REM) * VecHelpH(:NVEC_REM)
400 : ELSE ! real
401 : IF (set0) THEN
402 : smat%data_c(kj_off:kj_vec,ikG0) = chi*cfac(:NVEC_REM) * VecHelpS(:NVEC_REM)
403 : ELSE
404 : smat%data_c(kj_off:kj_vec,ikG0) = &
405 : smat%data_c(kj_off:kj_vec,ikG0) + chi*cfac(:NVEC_REM) * VecHelpS(:NVEC_REM)
406 : END IF
407 : hmat%data_c(kj_off:kj_vec,ikG0) = &
408 : hmat%data_c(kj_off:kj_vec,ikG0) + chi*cfac(:NVEC_REM) * VecHelpH(:NVEC_REM)
409 : END IF ! real
410 : END DO ! kj_off
411 : END DO ! ikG
412 : !$OMP END DO
413 : DEALLOCATE(plegend)
414 : DEALLOCATE(VecHelpS,VecHelpH)
415 : DEALLOCATE(cph_re,cph_im,cfac)
416 : DEALLOCATE(dot,fct,fct2)
417 : DEALLOCATE(xlegend)
418 : !$OMP END PARALLEL
419 17504 : CALL timestop("spherical setup")
420 17504 : RETURN
421 17504 : END SUBROUTINE hsmt_sph_cpu
422 :
423 : END MODULE m_hsmt_sph
|
