(git:15c7cac)
Loading...
Searching...
No Matches
construct_shg.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Calculation of the integrals over solid harmonic Gaussian(SHG) functions.
10!> Routines for (a|O(r12)|b) and overlap integrals (ab), (aba) and (abb).
11!> \par Literature (partly)
12!> T.J. Giese and D. M. York, J. Chem. Phys, 128, 064104 (2008)
13!> T. Helgaker, P Joergensen, J. Olsen, Molecular Electronic-Structure
14!> Theory, Wiley
15!> \par History
16!> created [04.2015]
17!> \author Dorothea Golze
18! **************************************************************************************************
19MODULE construct_shg
20 USE kinds, ONLY: dp
21 USE mathconstants, ONLY: dfac,&
22 fac
23 USE orbital_pointers, ONLY: indso,&
24 indso_inv,&
25 nsoset
26#include "../base/base_uses.f90"
27
28 IMPLICIT NONE
29
30 PRIVATE
31
32 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'construct_shg'
33
34! *** Public subroutines ***
35 PUBLIC :: get_real_scaled_solid_harmonic, get_w_matrix, get_dw_matrix, &
36 construct_int_shg_ab, construct_dev_shg_ab, construct_overlap_shg_aba, &
37 dev_overlap_shg_aba, construct_overlap_shg_abb, dev_overlap_shg_abb
38
39CONTAINS
40
41! **************************************************************************************************
42!> \brief computes the real scaled solid harmonics Rlm up to a given l
43!> \param Rlm_c cosine part of real scaled soldi harmonics
44!> \param Rlm_s sine part of real scaled soldi harmonics
45!> \param l maximal l quantum up to where Rlm is calculated
46!> \param r distance vector between a and b
47!> \param r2 square of distance vector
48! **************************************************************************************************
49 SUBROUTINE get_real_scaled_solid_harmonic(Rlm_c, Rlm_s, l, r, r2)
50
51 INTEGER, INTENT(IN) :: l
52 REAL(kind=dp), DIMENSION(0:l, -2*l:2*l), &
53 INTENT(OUT) :: rlm_s, rlm_c
54 REAL(kind=dp), DIMENSION(3) :: r
55 REAL(kind=dp) :: r2
56
57 INTEGER :: li, mi, prefac
58 REAL(kind=dp) :: rc, rc_00, rlm, rmlm, rplm, rs, rs_00, &
59 temp_c
60
61 rc_00 = 1.0_dp
62 rs_00 = 0.0_dp
63
64 rlm_c(0, 0) = rc_00
65 rlm_s(0, 0) = rs_00
66
67 ! generate elements Rmm
68 ! start
69 IF (l > 0) THEN
70 rc = -0.5_dp*r(1)*rc_00
71 rs = -0.5_dp*r(2)*rc_00
72 rlm_c(1, 1) = rc
73 rlm_s(1, 1) = rs
74 rlm_c(1, -1) = -rc
75 rlm_s(1, -1) = rs
76 END IF
77 DO li = 2, l
78 temp_c = (-r(1)*rc + r(2)*rs)/(real(2*(li - 1) + 2, dp))
79 rs = (-r(2)*rc - r(1)*rs)/(real(2*(li - 1) + 2, dp))
80 rc = temp_c
81 rlm_c(li, li) = rc
82 rlm_s(li, li) = rs
83 IF (modulo(li, 2) /= 0) THEN
84 rlm_c(li, -li) = -rc
85 rlm_s(li, -li) = rs
86 ELSE
87 rlm_c(li, -li) = rc
88 rlm_s(li, -li) = -rs
89 END IF
90 END DO
91
92 DO mi = 0, l - 1
93 rmlm = rlm_c(mi, mi)
94 rlm = r(3)*rlm_c(mi, mi)
95 rlm_c(mi + 1, mi) = rlm
96 IF (modulo(mi, 2) /= 0) THEN
97 rlm_c(mi + 1, -mi) = -rlm
98 ELSE
99 rlm_c(mi + 1, -mi) = rlm
100 END IF
101 DO li = mi + 2, l
102 prefac = (li + mi)*(li - mi)
103 rplm = (real(2*li - 1, dp)*r(3)*rlm - r2*rmlm)/real(prefac, dp)
104 rmlm = rlm
105 rlm = rplm
106 rlm_c(li, mi) = rlm
107 IF (modulo(mi, 2) /= 0) THEN
108 rlm_c(li, -mi) = -rlm
109 ELSE
110 rlm_c(li, -mi) = rlm
111 END IF
112 END DO
113 END DO
114 DO mi = 1, l - 1
115 rmlm = rlm_s(mi, mi)
116 rlm = r(3)*rlm_s(mi, mi)
117 rlm_s(mi + 1, mi) = rlm
118 IF (modulo(mi, 2) /= 0) THEN
119 rlm_s(mi + 1, -mi) = rlm
120 ELSE
121 rlm_s(mi + 1, -mi) = -rlm
122 END IF
123 DO li = mi + 2, l
124 prefac = (li + mi)*(li - mi)
125 rplm = (real(2*li - 1, dp)*r(3)*rlm - r2*rmlm)/real(prefac, dp)
126 rmlm = rlm
127 rlm = rplm
128 rlm_s(li, mi) = rlm
129 IF (modulo(mi, 2) /= 0) THEN
130 rlm_s(li, -mi) = rlm
131 ELSE
132 rlm_s(li, -mi) = -rlm
133 END IF
134 END DO
135 END DO
136
137 END SUBROUTINE get_real_scaled_solid_harmonic
138
139! **************************************************************************************************
140!> \brief Calculate the prefactor A(l,m) = (-1)^m \sqrt[(2-delta(m,0))(l+m)!(l-m)!]
141!> \param lmax maximal l quantum number
142!> \param A matrix storing the prefactor for a given l and m
143! **************************************************************************************************
144 SUBROUTINE get_alm(lmax, A)
145
146 INTEGER, INTENT(IN) :: lmax
147 REAL(kind=dp), DIMENSION(0:lmax, 0:lmax), &
148 INTENT(INOUT) :: a
149
150 INTEGER :: l, m
151 REAL(kind=dp) :: temp
152
153 DO l = 0, lmax
154 DO m = 0, l
155 temp = sqrt(fac(l + m)*fac(l - m))
156 IF (modulo(m, 2) /= 0) temp = -temp
157 IF (m /= 0) temp = temp*sqrt(2.0_dp)
158 a(l, m) = temp
159 END DO
160 END DO
161
162 END SUBROUTINE get_alm
163
164! **************************************************************************************************
165!> \brief calculates the prefactors for the derivatives of the W matrix
166!> \param lmax maximal l quantum number
167!> \param dA_p = A(l,m)/A(l-1,m+1)
168!> \param dA_m = A(l,m)/A(l-1,m-1)
169!> \param dA = A(l,m)/A(l-1,m)
170!> \note for m=0, W_l-1,-1 can't be read from Waux_mat, but we use
171!> W_l-1,-1 = -W_l-1,1 [cc(1), cs(2)] or W_l-1,-1 = W_l-1,1 [[sc(3), ss(4)], i.e.
172!> effectively we multiply dA_p by 2
173! **************************************************************************************************
174 SUBROUTINE get_da_prefactors(lmax, dA_p, dA_m, dA)
175
176 INTEGER, INTENT(IN) :: lmax
177 REAL(kind=dp), DIMENSION(0:lmax, 0:lmax), &
178 INTENT(INOUT) :: da_p, da_m, da
179
180 INTEGER :: l, m
181 REAL(kind=dp) :: bm, bm_m, bm_p
182
183 DO l = 0, lmax
184 DO m = 0, l
185 bm = 1.0_dp
186 bm_m = 1.0_dp
187 bm_p = 1.0_dp
188 IF (m /= 0) bm = sqrt(2.0_dp)
189 IF (m - 1 /= 0) bm_m = sqrt(2.0_dp)
190 IF (m + 1 /= 0) bm_p = sqrt(2.0_dp)
191 da_p(l, m) = -bm/bm_p*sqrt(real((l - m)*(l - m - 1), dp))
192 da_m(l, m) = -bm/bm_m*sqrt(real((l + m)*(l + m - 1), dp))
193 da(l, m) = 2.0_dp*sqrt(real((l + m)*(l - m), dp))
194 IF (m == 0) da_p(l, m) = 2.0_dp*da_p(l, m)
195 END DO
196 END DO
197 END SUBROUTINE get_da_prefactors
198
199! **************************************************************************************************
200!> \brief calculates the angular dependent-part of the SHG integrals,
201!> transformation matrix W, see literature above
202!> \param lamax array of maximal l quantum number on a;
203!> lamax(lb) with lb= 0..lbmax
204!> \param lbmax maximal l quantum number on b
205!> \param lmax maximal l quantum number
206!> \param Rc cosine part of real scaled solid harmonics
207!> \param Rs sine part of real scaled solid harmonics
208!> \param Waux_mat stores the angular-dependent part of the SHG integrals
209! **************************************************************************************************
210 SUBROUTINE get_w_matrix(lamax, lbmax, lmax, Rc, Rs, Waux_mat)
211
212 INTEGER, DIMENSION(:), POINTER :: lamax
213 INTEGER, INTENT(IN) :: lbmax, lmax
214 REAL(kind=dp), DIMENSION(0:lmax, -2*lmax:2*lmax), &
215 INTENT(IN) :: rc, rs
216 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: waux_mat
217
218 INTEGER :: j, k, la, labmin, laj, lb, lbj, ma, &
219 ma_m, ma_p, mb, mb_m, mb_p, nla, nlb
220 REAL(kind=dp) :: a_jk, a_lama, a_lbmb, alm_fac, delta_k, prefac, rca_m, rca_p, rcb_m, rcb_p, &
221 rsa_m, rsa_p, rsb_m, rsb_p, sign_fac, wa(4), wb(4), wmat(4)
222 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: a
223
224 wa(:) = 0.0_dp
225 wb(:) = 0.0_dp
226 wmat(:) = 0.0_dp
227
228 ALLOCATE (a(0:lmax, 0:lmax))
229 CALL get_alm(lmax, a)
230
231 DO lb = 0, lbmax
232 nlb = nsoset(lb - 1)
233 DO la = 0, lamax(lb)
234 nla = nsoset(la - 1)
235 labmin = min(la, lb)
236 DO mb = 0, lb
237 a_lbmb = a(lb, mb)
238 IF (modulo(lb, 2) /= 0) a_lbmb = -a_lbmb
239 DO ma = 0, la
240 a_lama = a(la, ma)
241 alm_fac = a_lama*a_lbmb
242 DO j = 0, labmin
243 laj = la - j
244 lbj = lb - j
245 prefac = alm_fac*real(2**(la + lb - j), dp)*dfac(2*j - 1)
246 delta_k = 0.5_dp
247 wmat = 0.0_dp
248 DO k = 0, j
249 ma_m = ma - k
250 ma_p = ma + k
251 IF (laj < abs(ma_m) .AND. laj < abs(ma_p)) cycle
252 mb_m = mb - k
253 mb_p = mb + k
254 IF (lbj < abs(mb_m) .AND. lbj < abs(mb_p)) cycle
255 IF (k /= 0) delta_k = 1.0_dp
256 a_jk = fac(j + k)*fac(j - k)
257 IF (k /= 0) a_jk = 2.0_dp*a_jk
258 IF (modulo(k, 2) /= 0) THEN
259 sign_fac = -1.0_dp
260 ELSE
261 sign_fac = 1.0_dp
262 END IF
263 rca_m = rc(laj, ma_m)
264 rsa_m = rs(laj, ma_m)
265 rca_p = rc(laj, ma_p)
266 rsa_p = rs(laj, ma_p)
267 rcb_m = rc(lbj, mb_m)
268 rsb_m = rs(lbj, mb_m)
269 rcb_p = rc(lbj, mb_p)
270 rsb_p = rs(lbj, mb_p)
271 wa(1) = delta_k*(rca_m + sign_fac*rca_p)
272 wb(1) = delta_k*(rcb_m + sign_fac*rcb_p)
273 wa(2) = -rsa_m + sign_fac*rsa_p
274 wb(2) = -rsb_m + sign_fac*rsb_p
275 wmat(1) = wmat(1) + prefac/a_jk*(wa(1)*wb(1) + wa(2)*wb(2))
276 IF (mb > 0) THEN
277 wb(3) = delta_k*(rsb_m + sign_fac*rsb_p)
278 wb(4) = rcb_m - sign_fac*rcb_p
279 wmat(2) = wmat(2) + prefac/a_jk*(wa(1)*wb(3) + wa(2)*wb(4))
280 END IF
281 IF (ma > 0) THEN
282 wa(3) = delta_k*(rsa_m + sign_fac*rsa_p)
283 wa(4) = rca_m - sign_fac*rca_p
284 wmat(3) = wmat(3) + prefac/a_jk*(wa(3)*wb(1) + wa(4)*wb(2))
285 END IF
286 IF (ma > 0 .AND. mb > 0) THEN
287 wmat(4) = wmat(4) + prefac/a_jk*(wa(3)*wb(3) + wa(4)*wb(4))
288 END IF
289 END DO
290 waux_mat(j + 1, nla + la + 1 + ma, nlb + lb + 1 + mb) = wmat(1)
291 IF (mb > 0) waux_mat(j + 1, nla + la + 1 + ma, nlb + lb + 1 - mb) = wmat(2)
292 IF (ma > 0) waux_mat(j + 1, nla + la + 1 - ma, nlb + lb + 1 + mb) = wmat(3)
293 IF (ma > 0 .AND. mb > 0) waux_mat(j + 1, nla + la + 1 - ma, nlb + lb + 1 - mb) = wmat(4)
294 END DO
295 END DO
296 END DO
297 END DO
298 END DO
299
300 DEALLOCATE (a)
301
302 END SUBROUTINE get_w_matrix
303
304! **************************************************************************************************
305!> \brief calculates derivatives of transformation matrix W,
306!> \param lamax array of maximal l quantum number on a;
307!> lamax(lb) with lb= 0..lbmax
308!> \param lbmax maximal l quantum number on b
309!> \param Waux_mat stores the angular-dependent part of the SHG integrals
310!> \param dWaux_mat stores the derivatives of the angular-dependent part of
311!> the SHG integrals
312! **************************************************************************************************
313 SUBROUTINE get_dw_matrix(lamax, lbmax, Waux_mat, dWaux_mat)
314
315 INTEGER, DIMENSION(:), POINTER :: lamax
316 INTEGER, INTENT(IN) :: lbmax
317 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: waux_mat
318 REAL(kind=dp), DIMENSION(:, :, :, :), &
319 INTENT(INOUT) :: dwaux_mat
320
321 INTEGER :: ima, imam, imb, imbm, ipa, ipam, ipb, &
322 ipbm, j, jmax, la, labm, labmin, lamb, &
323 lb, lmax, ma, mb, nla, nlam, nlb, nlbm
324 REAL(kind=dp) :: daa, daa_m, daa_p, dab, dab_m, dab_p
325 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: da, da_m, da_p, wam, wamm, wamp, wbm, &
326 wbmm, wbmp
327
328 jmax = min(maxval(lamax), lbmax)
329 ALLOCATE (wam(0:jmax, 4), wamm(0:jmax, 4), wamp(0:jmax, 4))
330 ALLOCATE (wbm(0:jmax, 4), wbmm(0:jmax, 4), wbmp(0:jmax, 4))
331
332 !*** get dA_p=A(l,m)/A(l-1,m+1)
333 !*** get dA_m=A(l,m)/A(l-1,m-1)
334 !*** get dA=2*A(l,m)/A(l-1,m)
335 lmax = max(maxval(lamax), lbmax)
336 ALLOCATE (da_p(0:lmax, 0:lmax), da_m(0:lmax, 0:lmax), da(0:lmax, 0:lmax))
337 CALL get_da_prefactors(lmax, da_p, da_m, da)
338
339 DO lb = 0, lbmax
340 nlb = nsoset(lb - 1)
341 nlbm = 0
342 IF (lb > 0) nlbm = nsoset(lb - 2)
343 DO la = 0, lamax(lb)
344 nla = nsoset(la - 1)
345 nlam = 0
346 IF (la > 0) nlam = nsoset(la - 2)
347 labmin = min(la, lb)
348 lamb = min(la - 1, lb)
349 labm = min(la, lb - 1)
350 DO mb = 0, lb
351 dab = da(lb, mb)
352 dab_p = da_p(lb, mb)
353 dab_m = da_m(lb, mb)
354 ipb = nlb + lb + mb + 1
355 imb = nlb + lb - mb + 1
356 ipbm = nlbm + lb + mb
357 imbm = nlbm + lb - mb
358 DO ma = 0, la
359 daa = da(la, ma)
360 daa_p = da_p(la, ma)
361 daa_m = da_m(la, ma)
362 ipa = nla + la + ma + 1
363 ima = nla + la - ma + 1
364 ipam = nlam + la + ma
365 imam = nlam + la - ma
366 wam(:, :) = 0.0_dp
367 wamm(:, :) = 0.0_dp
368 wamp(:, :) = 0.0_dp
369 !*** Wam: la-1, ma
370 IF (ma <= la - 1) THEN
371 wam(0:lamb, 1) = waux_mat(1:lamb + 1, ipam, ipb)
372 IF (mb > 0) wam(0:lamb, 2) = waux_mat(1:lamb + 1, ipam, imb)
373 IF (ma > 0) wam(0:lamb, 3) = waux_mat(1:lamb + 1, imam, ipb)
374 IF (ma > 0 .AND. mb > 0) wam(0:lamb, 4) = waux_mat(1:lamb + 1, imam, imb)
375 END IF
376 !*** Wamm: la-1, ma-1
377 IF (ma - 1 >= 0) THEN
378 wamm(0:lamb, 1) = waux_mat(1:lamb + 1, ipam - 1, ipb)
379 IF (mb > 0) wamm(0:lamb, 2) = waux_mat(1:lamb + 1, ipam - 1, imb)
380 ! order: e.g. -1 0 1, if < 0 |m|, -1 means -m+1
381 IF (ma - 1 > 0) wamm(0:lamb, 3) = waux_mat(1:lamb + 1, imam + 1, ipb)
382 IF (ma - 1 > 0 .AND. mb > 0) wamm(0:lamb, 4) = waux_mat(1:lamb + 1, imam + 1, imb)
383 END IF
384 !*** Wamp: la-1, ma+1
385 IF (ma + 1 <= la - 1) THEN
386 wamp(0:lamb, 1) = waux_mat(1:lamb + 1, ipam + 1, ipb)
387 IF (mb > 0) wamp(0:lamb, 2) = waux_mat(1:lamb + 1, ipam + 1, imb)
388 IF (ma + 1 > 0) wamp(0:lamb, 3) = waux_mat(1:lamb + 1, imam - 1, ipb)
389 IF (ma + 1 > 0 .AND. mb > 0) wamp(0:lamb, 4) = waux_mat(1:lamb + 1, imam - 1, imb)
390 END IF
391 wbm(:, :) = 0.0_dp
392 wbmm(:, :) = 0.0_dp
393 wbmp(:, :) = 0.0_dp
394 !*** Wbm: lb-1, mb
395 IF (mb <= lb - 1) THEN
396 wbm(0:labm, 1) = waux_mat(1:labm + 1, ipa, ipbm)
397 IF (mb > 0) wbm(0:labm, 2) = waux_mat(1:labm + 1, ipa, imbm)
398 IF (ma > 0) wbm(0:labm, 3) = waux_mat(1:labm + 1, ima, ipbm)
399 IF (ma > 0 .AND. mb > 0) wbm(0:labm, 4) = waux_mat(1:labm + 1, ima, imbm)
400 END IF
401 !*** Wbmm: lb-1, mb-1
402 IF (mb - 1 >= 0) THEN
403 wbmm(0:labm, 1) = waux_mat(1:labm + 1, ipa, ipbm - 1)
404 IF (mb - 1 > 0) wbmm(0:labm, 2) = waux_mat(1:labm + 1, ipa, imbm + 1)
405 IF (ma > 0) wbmm(0:labm, 3) = waux_mat(1:labm + 1, ima, ipbm - 1)
406 IF (ma > 0 .AND. mb - 1 > 0) wbmm(0:labm, 4) = waux_mat(1:labm + 1, ima, imbm + 1)
407 END IF
408 !*** Wbmp: lb-1, mb+1
409 IF (mb + 1 <= lb - 1) THEN
410 wbmp(0:labm, 1) = waux_mat(1:labm + 1, ipa, ipbm + 1)
411 IF (mb + 1 > 0) wbmp(0:labm, 2) = waux_mat(1:labm + 1, ipa, imbm - 1)
412 IF (ma > 0) wbmp(0:labm, 3) = waux_mat(1:labm + 1, ima, ipbm + 1)
413 IF (ma > 0 .AND. mb + 1 > 0) wbmp(0:labm, 4) = waux_mat(1:labm + 1, ima, imbm - 1)
414 END IF
415 DO j = 0, labmin
416 !*** x component
417 dwaux_mat(1, j + 1, ipa, ipb) = daa_p*wamp(j, 1) - daa_m*wamm(j, 1) &
418 - dab_p*wbmp(j, 1) + dab_m*wbmm(j, 1)
419 IF (mb > 0) THEN
420 dwaux_mat(1, j + 1, ipa, imb) = daa_p*wamp(j, 2) - daa_m*wamm(j, 2) &
421 - dab_p*wbmp(j, 2) + dab_m*wbmm(j, 2)
422 END IF
423 IF (ma > 0) THEN
424 dwaux_mat(1, j + 1, ima, ipb) = daa_p*wamp(j, 3) - daa_m*wamm(j, 3) &
425 - dab_p*wbmp(j, 3) + dab_m*wbmm(j, 3)
426 END IF
427 IF (ma > 0 .AND. mb > 0) THEN
428 dwaux_mat(1, j + 1, ima, imb) = daa_p*wamp(j, 4) - daa_m*wamm(j, 4) &
429 - dab_p*wbmp(j, 4) + dab_m*wbmm(j, 4)
430 END IF
431
432 !**** y component
433 dwaux_mat(2, j + 1, ipa, ipb) = daa_p*wamp(j, 3) + daa_m*wamm(j, 3) &
434 - dab_p*wbmp(j, 2) - dab_m*wbmm(j, 2)
435 IF (mb > 0) THEN
436 dwaux_mat(2, j + 1, ipa, imb) = daa_p*wamp(j, 4) + daa_m*wamm(j, 4) &
437 + dab_p*wbmp(j, 1) + dab_m*wbmm(j, 1)
438 END IF
439 IF (ma > 0) THEN
440 dwaux_mat(2, j + 1, ima, ipb) = -daa_p*wamp(j, 1) - daa_m*wamm(j, 1) &
441 - dab_p*wbmp(j, 4) - dab_m*wbmm(j, 4)
442 END IF
443 IF (ma > 0 .AND. mb > 0) THEN
444 dwaux_mat(2, j + 1, ima, imb) = -daa_p*wamp(j, 2) - daa_m*wamm(j, 2) &
445 + dab_p*wbmp(j, 3) + dab_m*wbmm(j, 3)
446 END IF
447 !**** z compnent
448 dwaux_mat(3, j + 1, ipa, ipb) = daa*wam(j, 1) - dab*wbm(j, 1)
449 IF (mb > 0) THEN
450 dwaux_mat(3, j + 1, ipa, imb) = daa*wam(j, 2) - dab*wbm(j, 2)
451 END IF
452 IF (ma > 0) THEN
453 dwaux_mat(3, j + 1, ima, ipb) = daa*wam(j, 3) - dab*wbm(j, 3)
454 END IF
455 IF (ma > 0 .AND. mb > 0) THEN
456 dwaux_mat(3, j + 1, ima, imb) = daa*wam(j, 4) - dab*wbm(j, 4)
457 END IF
458
459 END DO
460 END DO
461 END DO
462 END DO
463 END DO
464
465 DEALLOCATE (wam, wamm, wamp)
466 DEALLOCATE (wbm, wbmm, wbmp)
467 DEALLOCATE (da, da_p, da_m)
468
469 END SUBROUTINE get_dw_matrix
470
471! **************************************************************************************************
472!> \brief calculates [ab] SHG overlap integrals using precomputed angular-
473!> dependent part
474!> \param la set of l quantum number on a
475!> \param first_sgfa indexing
476!> \param nshella number of shells for a
477!> \param lb set of l quantum number on b
478!> \param first_sgfb indexing
479!> \param nshellb number of shells for b
480!> \param swork_cont contracted and normalized [s|s] integrals
481!> \param Waux_mat precomputed angular-dependent part
482!> \param sab contracted integral of spherical harmonic Gaussianslm
483! **************************************************************************************************
484 SUBROUTINE construct_int_shg_ab(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
485 swork_cont, Waux_mat, sab)
486
487 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
488 INTEGER, INTENT(IN) :: nshella
489 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
490 INTEGER, INTENT(IN) :: nshellb
491 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: swork_cont, waux_mat
492 REAL(kind=dp), DIMENSION(:, :), INTENT(INOUT) :: sab
493
494 INTEGER :: fnla, fnlb, fsgfa, fsgfb, ishella, j, &
495 jshellb, labmin, lai, lbj, lnla, lnlb, &
496 lsgfa, lsgfb, mai, mbj
497 REAL(kind=dp) :: prefac
498
499 DO jshellb = 1, nshellb
500 lbj = lb(jshellb)
501 fnlb = nsoset(lbj - 1) + 1
502 lnlb = nsoset(lbj)
503 fsgfb = first_sgfb(jshellb)
504 lsgfb = fsgfb + 2*lbj
505 DO ishella = 1, nshella
506 lai = la(ishella)
507 fnla = nsoset(lai - 1) + 1
508 lnla = nsoset(lai)
509 fsgfa = first_sgfa(ishella)
510 lsgfa = fsgfa + 2*lai
511 labmin = min(lai, lbj)
512 DO mbj = 0, 2*lbj
513 DO mai = 0, 2*lai
514 DO j = 0, labmin
515 prefac = swork_cont(lai + lbj - j + 1, ishella, jshellb)
516 sab(fsgfa + mai, fsgfb + mbj) = sab(fsgfa + mai, fsgfb + mbj) &
517 + prefac*waux_mat(j + 1, fnla + mai, fnlb + mbj)
518 END DO
519 END DO
520 END DO
521 END DO
522 END DO
523
524 END SUBROUTINE construct_int_shg_ab
525
526! **************************************************************************************************
527!> \brief calculates derivatives of [ab] SHG overlap integrals using precomputed
528!> angular-dependent part
529!> \param la set of l quantum number on a
530!> \param first_sgfa indexing
531!> \param nshella number of shells for a
532!> \param lb set of l quantum number on b
533!> \param first_sgfb indexing
534!> \param nshellb number of shells for b
535!> \param rab distance vector Ra-Rb
536!> \param swork_cont contracted and normalized [s|s] integrals
537!> \param Waux_mat precomputed angular-dependent part
538!> \param dWaux_mat ...
539!> \param dsab derivative of contracted integral of spherical harmonic Gaussians
540! **************************************************************************************************
541 SUBROUTINE construct_dev_shg_ab(la, first_sgfa, nshella, lb, first_sgfb, nshellb, rab, &
542 swork_cont, Waux_mat, dWaux_mat, dsab)
543
544 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
545 INTEGER, INTENT(IN) :: nshella
546 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
547 INTEGER, INTENT(IN) :: nshellb
548 REAL(kind=dp), INTENT(IN) :: rab(3)
549 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: swork_cont, waux_mat
550 REAL(kind=dp), DIMENSION(:, :, :, :), INTENT(IN) :: dwaux_mat
551 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: dsab
552
553 INTEGER :: fnla, fnlb, fsgfa, fsgfb, i, ishella, j, &
554 jshellb, labmin, lai, lbj, lnla, lnlb, &
555 lsgfa, lsgfb
556 REAL(kind=dp) :: dprefac, prefac, rabx2(3)
557
558 rabx2(:) = 2.0_dp*rab
559 DO jshellb = 1, nshellb
560 lbj = lb(jshellb)
561 fnlb = nsoset(lbj - 1) + 1
562 lnlb = nsoset(lbj)
563 fsgfb = first_sgfb(jshellb)
564 lsgfb = fsgfb + 2*lbj
565 DO ishella = 1, nshella
566 lai = la(ishella)
567 fnla = nsoset(lai - 1) + 1
568 lnla = nsoset(lai)
569 fsgfa = first_sgfa(ishella)
570 lsgfa = fsgfa + 2*lai
571 labmin = min(lai, lbj)
572 DO j = 0, labmin
573 prefac = swork_cont(lai + lbj - j + 1, ishella, jshellb)
574 dprefac = swork_cont(lai + lbj - j + 2, ishella, jshellb) !j+1
575 DO i = 1, 3
576 dsab(fsgfa:lsgfa, fsgfb:lsgfb, i) = dsab(fsgfa:lsgfa, fsgfb:lsgfb, i) &
577 + rabx2(i)*dprefac*waux_mat(j + 1, fnla:lnla, fnlb:lnlb) &
578 + prefac*dwaux_mat(i, j + 1, fnla:lnla, fnlb:lnlb)
579 END DO
580 END DO
581 END DO
582 END DO
583
584 END SUBROUTINE construct_dev_shg_ab
585
586! **************************************************************************************************
587!> \brief calculates [aba] SHG overlap integrals using precomputed angular-
588!> dependent part
589!> \param la set of l quantum number on a, orbital basis
590!> \param first_sgfa indexing
591!> \param nshella number of shells for a, orbital basis
592!> \param lb set of l quantum number on b. orbital basis
593!> \param first_sgfb indexing
594!> \param nshellb number of shells for b, orbital basis
595!> \param lca of l quantum number on a, aux basis
596!> \param first_sgfca indexing
597!> \param nshellca number of shells for a, aux basis
598!> \param cg_coeff Clebsch-Gordon coefficients
599!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
600!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
601!> \param swork_cont contracted and normalized [s|ra^n|s] integrals
602!> \param Waux_mat precomputed angular-dependent part
603!> \param saba contracted overlap [aba] of spherical harmonic Gaussians
604! **************************************************************************************************
605 SUBROUTINE construct_overlap_shg_aba(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
606 lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, &
607 ncg_none0, swork_cont, Waux_mat, saba)
608
609 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
610 INTEGER, INTENT(IN) :: nshella
611 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
612 INTEGER, INTENT(IN) :: nshellb
613 INTEGER, DIMENSION(:), INTENT(IN) :: lca, first_sgfca
614 INTEGER, INTENT(IN) :: nshellca
615 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
616 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
617 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
618 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
619 INTENT(IN) :: swork_cont
620 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
621 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: saba
622
623 INTEGER :: ia, il, ilist, ishella, isoa1, isoa2, isoaa, j, jb, jshellb, ka, kshella, laa, &
624 labmin, lai, lak, lbj, maa, mai, mak, mbj, nla, nlb, sgfa, sgfb, sgfca
625 REAL(kind=dp) :: prefac, stemp
626
627 DO kshella = 1, nshellca
628 lak = lca(kshella)
629 sgfca = first_sgfca(kshella)
630 ka = sgfca + lak
631 DO jshellb = 1, nshellb
632 lbj = lb(jshellb)
633 nlb = nsoset(lbj - 1) + lbj + 1
634 sgfb = first_sgfb(jshellb)
635 jb = sgfb + lbj
636 DO ishella = 1, nshella
637 lai = la(ishella)
638 sgfa = first_sgfa(ishella)
639 ia = sgfa + lai
640 DO mai = -lai, lai, 1
641 DO mak = -lak, lak, 1
642 isoa1 = indso_inv(lai, mai)
643 isoa2 = indso_inv(lak, mak)
644 DO mbj = -lbj, lbj, 1
645 DO ilist = 1, ncg_none0(isoa1, isoa2)
646 isoaa = cg_none0_list(isoa1, isoa2, ilist)
647 laa = indso(1, isoaa)
648 maa = indso(2, isoaa)
649 nla = nsoset(laa - 1) + laa + 1
650 labmin = min(laa, lbj)
651 il = int((lai + lak - laa)/2)
652 stemp = 0.0_dp
653 DO j = 0, labmin
654 prefac = swork_cont(laa + lbj - j + 1, il, ishella, jshellb, kshella)
655 stemp = stemp + prefac*waux_mat(j + 1, nla + maa, nlb + mbj)
656 END DO
657 saba(ia + mai, jb + mbj, ka + mak) = saba(ia + mai, jb + mbj, ka + mak) + cg_coeff(isoa1, isoa2, isoaa)*stemp
658 END DO
659 END DO
660 END DO
661 END DO
662 END DO
663 END DO
664 END DO
665
666 END SUBROUTINE construct_overlap_shg_aba
667
668! **************************************************************************************************
669!> \brief calculates derivatives of [aba] SHG overlap integrals using
670!> precomputed angular-dependent part
671!> \param la set of l quantum number on a, orbital basis
672!> \param first_sgfa indexing
673!> \param nshella number of shells for a, orbital basis
674!> \param lb set of l quantum number on b. orbital basis
675!> \param first_sgfb indexing
676!> \param nshellb number of shells for b, orbital basis
677!> \param lca of l quantum number on a, aux basis
678!> \param first_sgfca indexing
679!> \param nshellca number of shells for a, aux basis
680!> \param cg_coeff Clebsch-Gordon coefficients
681!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
682!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
683!> \param rab distance vector Ra-Rb
684!> \param swork_cont contracted and normalized [s|ra^n|s] integrals
685!> \param Waux_mat precomputed angular-dependent part
686!> \param dWaux_mat derivatives of precomputed angular-dependent part
687!> \param dsaba derivative of contracted overlap [aba] of spherical harmonic
688!> Gaussians
689! **************************************************************************************************
690 SUBROUTINE dev_overlap_shg_aba(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
691 lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, &
692 ncg_none0, rab, swork_cont, Waux_mat, dWaux_mat, dsaba)
693
694 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
695 INTEGER, INTENT(IN) :: nshella
696 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
697 INTEGER, INTENT(IN) :: nshellb
698 INTEGER, DIMENSION(:), INTENT(IN) :: lca, first_sgfca
699 INTEGER, INTENT(IN) :: nshellca
700 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
701 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
702 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
703 REAL(kind=dp), INTENT(IN) :: rab(3)
704 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
705 INTENT(IN) :: swork_cont
706 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
707 REAL(kind=dp), DIMENSION(:, :, :, :), INTENT(IN) :: dwaux_mat
708 REAL(kind=dp), DIMENSION(:, :, :, :), &
709 INTENT(INOUT) :: dsaba
710
711 INTEGER :: i, ia, il, ilist, ishella, isoa1, isoa2, isoaa, j, jb, jshellb, ka, kshella, laa, &
712 labmin, lai, lak, lbj, maa, mai, mak, mbj, nla, nlb, sgfa, sgfb, sgfca
713 REAL(kind=dp) :: dprefac, dtemp(3), prefac, rabx2(3)
714
715 rabx2(:) = 2.0_dp*rab
716
717 DO kshella = 1, nshellca
718 lak = lca(kshella)
719 sgfca = first_sgfca(kshella)
720 ka = sgfca + lak
721 DO jshellb = 1, nshellb
722 lbj = lb(jshellb)
723 nlb = nsoset(lbj - 1) + lbj + 1
724 sgfb = first_sgfb(jshellb)
725 jb = sgfb + lbj
726 DO ishella = 1, nshella
727 lai = la(ishella)
728 sgfa = first_sgfa(ishella)
729 ia = sgfa + lai
730 DO mai = -lai, lai, 1
731 DO mak = -lak, lak, 1
732 isoa1 = indso_inv(lai, mai)
733 isoa2 = indso_inv(lak, mak)
734 DO mbj = -lbj, lbj, 1
735 DO ilist = 1, ncg_none0(isoa1, isoa2)
736 isoaa = cg_none0_list(isoa1, isoa2, ilist)
737 laa = indso(1, isoaa)
738 maa = indso(2, isoaa)
739 nla = nsoset(laa - 1) + laa + 1
740 labmin = min(laa, lbj)
741 il = (lai + lak - laa)/2 ! lai+lak-laa always even
742 dtemp = 0.0_dp
743 DO j = 0, labmin
744 prefac = swork_cont(laa + lbj - j + 1, il, ishella, jshellb, kshella)
745 dprefac = swork_cont(laa + lbj - j + 2, il, ishella, jshellb, kshella)
746 DO i = 1, 3
747 dtemp(i) = dtemp(i) + rabx2(i)*dprefac*waux_mat(j + 1, nla + maa, nlb + mbj) &
748 + prefac*dwaux_mat(i, j + 1, nla + maa, nlb + mbj)
749 END DO
750 END DO
751 DO i = 1, 3
752 dsaba(ia + mai, jb + mbj, ka + mak, i) = dsaba(ia + mai, jb + mbj, ka + mak, i) &
753 + cg_coeff(isoa1, isoa2, isoaa)*dtemp(i)
754 END DO
755 END DO
756 END DO
757 END DO
758 END DO
759 END DO
760 END DO
761 END DO
762
763 END SUBROUTINE dev_overlap_shg_aba
764
765! **************************************************************************************************
766!> \brief calculates [abb] SHG overlap integrals using precomputed angular-
767!> dependent part
768!> \param la set of l quantum number on a, orbital basis
769!> \param first_sgfa indexing
770!> \param nshella number of shells for a, orbital basis
771!> \param lb set of l quantum number on b. orbital basis
772!> \param first_sgfb indexing
773!> \param nshellb number of shells for b, orbital basis
774!> \param lcb l quantum number on b, aux basis
775!> \param first_sgfcb indexing
776!> \param nshellcb number of shells for b, aux basis
777!> \param cg_coeff Clebsch-Gordon coefficients
778!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
779!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
780!> \param swork_cont contracted and normalized [s|rb^n|s] integrals
781!> \param Waux_mat precomputed angular-dependent part
782!> \param sabb contracted overlap [abb] of spherical harmonic Gaussians
783! **************************************************************************************************
784 SUBROUTINE construct_overlap_shg_abb(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
785 lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, &
786 ncg_none0, swork_cont, Waux_mat, sabb)
787
788 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
789 INTEGER, INTENT(IN) :: nshella
790 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
791 INTEGER, INTENT(IN) :: nshellb
792 INTEGER, DIMENSION(:), INTENT(IN) :: lcb, first_sgfcb
793 INTEGER, INTENT(IN) :: nshellcb
794 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
795 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
796 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
797 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
798 INTENT(IN) :: swork_cont
799 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
800 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: sabb
801
802 INTEGER :: ia, il, ilist, ishella, isob1, isob2, isobb, j, jb, jshellb, kb, kshellb, labmin, &
803 lai, lbb, lbj, lbk, mai, mbb, mbj, mbk, nla, nlb, sgfa, sgfb, sgfcb
804 REAL(kind=dp) :: prefac, stemp, tsign
805
806 DO kshellb = 1, nshellcb
807 lbk = lcb(kshellb)
808 sgfcb = first_sgfcb(kshellb)
809 kb = sgfcb + lbk
810 DO jshellb = 1, nshellb
811 lbj = lb(jshellb)
812 sgfb = first_sgfb(jshellb)
813 jb = sgfb + lbj
814 DO ishella = 1, nshella
815 lai = la(ishella)
816 nla = nsoset(lai - 1) + lai + 1
817 sgfa = first_sgfa(ishella)
818 ia = sgfa + lai
819 DO mbj = -lbj, lbj, 1
820 DO mbk = -lbk, lbk, 1
821 isob1 = indso_inv(lbj, mbj)
822 isob2 = indso_inv(lbk, mbk)
823 DO mai = -lai, lai, 1
824 DO ilist = 1, ncg_none0(isob1, isob2)
825 isobb = cg_none0_list(isob1, isob2, ilist)
826 lbb = indso(1, isobb)
827 mbb = indso(2, isobb)
828 nlb = nsoset(lbb - 1) + lbb + 1
829 ! tsgin: because we take the transpose of auxmat (calculated for (la,lb))
830 tsign = 1.0_dp
831 IF (modulo(lbb - lai, 2) /= 0) tsign = -1.0_dp
832 labmin = min(lai, lbb)
833 il = int((lbj + lbk - lbb)/2)
834 stemp = 0.0_dp
835 DO j = 0, labmin
836 prefac = swork_cont(lai + lbb - j + 1, il, ishella, jshellb, kshellb)
837 stemp = stemp + prefac*waux_mat(j + 1, nlb + mbb, nla + mai)
838 END DO
839 sabb(ia + mai, jb + mbj, kb + mbk) = sabb(ia + mai, jb + mbj, kb + mbk) + tsign*cg_coeff(isob1, isob2, isobb)*stemp
840 END DO
841 END DO
842 END DO
843 END DO
844 END DO
845 END DO
846 END DO
847
848 END SUBROUTINE construct_overlap_shg_abb
849
850! **************************************************************************************************
851!> \brief calculates derivatives of [abb] SHG overlap integrals using
852!> precomputed angular-dependent part
853!> \param la set of l quantum number on a, orbital basis
854!> \param first_sgfa indexing
855!> \param nshella number of shells for a, orbital basis
856!> \param lb set of l quantum number on b. orbital basis
857!> \param first_sgfb indexing
858!> \param nshellb number of shells for b, orbital basis
859!> \param lcb l quantum number on b, aux basis
860!> \param first_sgfcb indexing
861!> \param nshellcb number of shells for b, aux basis
862!> \param cg_coeff Clebsch-Gordon coefficients
863!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
864!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
865!> \param rab distance vector Ra-Rb
866!> \param swork_cont contracted and normalized [s|rb^n|s] integrals
867!> \param Waux_mat precomputed angular-dependent part
868!> \param dWaux_mat derivatives of precomputed angular-dependent part
869!> \param dsabb derivative of contracted overlap [abb] of spherical harmonic
870!> Gaussians
871! **************************************************************************************************
872 SUBROUTINE dev_overlap_shg_abb(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
873 lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, &
874 ncg_none0, rab, swork_cont, Waux_mat, dWaux_mat, dsabb)
875
876 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
877 INTEGER, INTENT(IN) :: nshella
878 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
879 INTEGER, INTENT(IN) :: nshellb
880 INTEGER, DIMENSION(:), INTENT(IN) :: lcb, first_sgfcb
881 INTEGER, INTENT(IN) :: nshellcb
882 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
883 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
884 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
885 REAL(kind=dp), INTENT(IN) :: rab(3)
886 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
887 INTENT(IN) :: swork_cont
888 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
889 REAL(kind=dp), DIMENSION(:, :, :, :), INTENT(IN) :: dwaux_mat
890 REAL(kind=dp), DIMENSION(:, :, :, :), &
891 INTENT(INOUT) :: dsabb
892
893 INTEGER :: i, ia, il, ilist, ishella, isob1, isob2, isobb, j, jb, jshellb, kb, kshellb, &
894 labmin, lai, lbb, lbj, lbk, mai, mbb, mbj, mbk, nla, nlb, sgfa, sgfb, sgfcb
895 REAL(kind=dp) :: dprefac, dtemp(3), prefac, rabx2(3), &
896 tsign
897
898 rabx2(:) = 2.0_dp*rab
899
900 DO kshellb = 1, nshellcb
901 lbk = lcb(kshellb)
902 sgfcb = first_sgfcb(kshellb)
903 kb = sgfcb + lbk
904 DO jshellb = 1, nshellb
905 lbj = lb(jshellb)
906 sgfb = first_sgfb(jshellb)
907 jb = sgfb + lbj
908 DO ishella = 1, nshella
909 lai = la(ishella)
910 nla = nsoset(lai - 1) + lai + 1
911 sgfa = first_sgfa(ishella)
912 ia = sgfa + lai
913 DO mbj = -lbj, lbj, 1
914 DO mbk = -lbk, lbk, 1
915 isob1 = indso_inv(lbj, mbj)
916 isob2 = indso_inv(lbk, mbk)
917 DO mai = -lai, lai, 1
918 DO ilist = 1, ncg_none0(isob1, isob2)
919 isobb = cg_none0_list(isob1, isob2, ilist)
920 lbb = indso(1, isobb)
921 mbb = indso(2, isobb)
922 nlb = nsoset(lbb - 1) + lbb + 1
923 ! tsgin: because we take the transpose of auxmat (calculated for (la,lb))
924 tsign = 1.0_dp
925 IF (modulo(lbb - lai, 2) /= 0) tsign = -1.0_dp
926 labmin = min(lai, lbb)
927 il = (lbj + lbk - lbb)/2
928 dtemp = 0.0_dp
929 DO j = 0, labmin
930 prefac = swork_cont(lai + lbb - j + 1, il, ishella, jshellb, kshellb)
931 dprefac = swork_cont(lai + lbb - j + 2, il, ishella, jshellb, kshellb)
932 DO i = 1, 3
933 dtemp(i) = dtemp(i) + rabx2(i)*dprefac*waux_mat(j + 1, nlb + mbb, nla + mai) &
934 + prefac*dwaux_mat(i, j + 1, nlb + mbb, nla + mai)
935 END DO
936 END DO
937 DO i = 1, 3
938 dsabb(ia + mai, jb + mbj, kb + mbk, i) = dsabb(ia + mai, jb + mbj, kb + mbk, i) &
939 + tsign*cg_coeff(isob1, isob2, isobb)*dtemp(i)
940 END DO
941 END DO
942 END DO
943 END DO
944 END DO
945 END DO
946 END DO
947 END DO
948
949 END SUBROUTINE dev_overlap_shg_abb
950
951END MODULE construct_shg
952
modulo
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Definition grid_common.h:125
construct_shg
Calculation of the integrals over solid harmonic Gaussian(SHG) functions. Routines for (a|O(r12)|b) a...
Definition construct_shg.F:19
construct_shg::dev_overlap_shg_aba
subroutine, public dev_overlap_shg_aba(la, first_sgfa, nshella, lb, first_sgfb, nshellb, lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, ncg_none0, rab, swork_cont, waux_mat, dwaux_mat, dsaba)
calculates derivatives of [aba] SHG overlap integrals using precomputed angular-dependent part
Definition construct_shg.F:693
construct_shg::get_dw_matrix
subroutine, public get_dw_matrix(lamax, lbmax, waux_mat, dwaux_mat)
calculates derivatives of transformation matrix W,
Definition construct_shg.F:314
construct_shg::get_w_matrix
subroutine, public get_w_matrix(lamax, lbmax, lmax, rc, rs, waux_mat)
calculates the angular dependent-part of the SHG integrals, transformation matrix W,...
Definition construct_shg.F:211
construct_shg::construct_dev_shg_ab
subroutine, public construct_dev_shg_ab(la, first_sgfa, nshella, lb, first_sgfb, nshellb, rab, swork_cont, waux_mat, dwaux_mat, dsab)
calculates derivatives of [ab] SHG overlap integrals using precomputed angular-dependent part
Definition construct_shg.F:543
construct_shg::construct_overlap_shg_abb
subroutine, public construct_overlap_shg_abb(la, first_sgfa, nshella, lb, first_sgfb, nshellb, lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, ncg_none0, swork_cont, waux_mat, sabb)
calculates [abb] SHG overlap integrals using precomputed angular- dependent part
Definition construct_shg.F:787
construct_shg::construct_int_shg_ab
subroutine, public construct_int_shg_ab(la, first_sgfa, nshella, lb, first_sgfb, nshellb, swork_cont, waux_mat, sab)
calculates [ab] SHG overlap integrals using precomputed angular- dependent part
Definition construct_shg.F:486
construct_shg::dev_overlap_shg_abb
subroutine, public dev_overlap_shg_abb(la, first_sgfa, nshella, lb, first_sgfb, nshellb, lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, ncg_none0, rab, swork_cont, waux_mat, dwaux_mat, dsabb)
calculates derivatives of [abb] SHG overlap integrals using precomputed angular-dependent part
Definition construct_shg.F:875
construct_shg::construct_overlap_shg_aba
subroutine, public construct_overlap_shg_aba(la, first_sgfa, nshella, lb, first_sgfb, nshellb, lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, ncg_none0, swork_cont, waux_mat, saba)
calculates [aba] SHG overlap integrals using precomputed angular- dependent part
Definition construct_shg.F:608
construct_shg::get_real_scaled_solid_harmonic
subroutine, public get_real_scaled_solid_harmonic(rlm_c, rlm_s, l, r, r2)
computes the real scaled solid harmonics Rlm up to a given l
Definition construct_shg.F:50
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
mathconstants
Definition of mathematical constants and functions.
Definition mathconstants.F:16
mathconstants::dfac
real(kind=dp), dimension(-1:2 *maxfac+1), parameter, public dfac
Definition mathconstants.F:61
mathconstants::fac
real(kind=dp), dimension(0:maxfac), parameter, public fac
Definition mathconstants.F:37
orbital_pointers
Provides Cartesian and spherical orbital pointers and indices.
Definition orbital_pointers.F:15
orbital_pointers::nsoset
integer, dimension(:), allocatable, public nsoset
Definition orbital_pointers.F:42
orbital_pointers::indso
integer, dimension(:, :), allocatable, public indso
Definition orbital_pointers.F:43
orbital_pointers::indso_inv
integer, dimension(:, :), allocatable, public indso_inv
Definition orbital_pointers.F:43