(git:374b731)
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-2024 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 IF (ma - 1 > 0) wamm(0:lamb, 3) = waux_mat(1:lamb + 1, imam + 1, ipb) !order: e.g. -1 0 1, if < 0 |m|, -1 means -m+1
381 IF (ma - 1 > 0 .AND. mb > 0) wamm(0:lamb, 4) = waux_mat(1:lamb + 1, imam + 1, imb)
382 END IF
383 !*** Wamp: la-1, ma+1
384 IF (ma + 1 <= la - 1) THEN
385 wamp(0:lamb, 1) = waux_mat(1:lamb + 1, ipam + 1, ipb)
386 IF (mb > 0) wamp(0:lamb, 2) = waux_mat(1:lamb + 1, ipam + 1, imb)
387 IF (ma + 1 > 0) wamp(0:lamb, 3) = waux_mat(1:lamb + 1, imam - 1, ipb)
388 IF (ma + 1 > 0 .AND. mb > 0) wamp(0:lamb, 4) = waux_mat(1:lamb + 1, imam - 1, imb)
389 END IF
390 wbm(:, :) = 0.0_dp
391 wbmm(:, :) = 0.0_dp
392 wbmp(:, :) = 0.0_dp
393 !*** Wbm: lb-1, mb
394 IF (mb <= lb - 1) THEN
395 wbm(0:labm, 1) = waux_mat(1:labm + 1, ipa, ipbm)
396 IF (mb > 0) wbm(0:labm, 2) = waux_mat(1:labm + 1, ipa, imbm)
397 IF (ma > 0) wbm(0:labm, 3) = waux_mat(1:labm + 1, ima, ipbm)
398 IF (ma > 0 .AND. mb > 0) wbm(0:labm, 4) = waux_mat(1:labm + 1, ima, imbm)
399 END IF
400 !*** Wbmm: lb-1, mb-1
401 IF (mb - 1 >= 0) THEN
402 wbmm(0:labm, 1) = waux_mat(1:labm + 1, ipa, ipbm - 1)
403 IF (mb - 1 > 0) wbmm(0:labm, 2) = waux_mat(1:labm + 1, ipa, imbm + 1)
404 IF (ma > 0) wbmm(0:labm, 3) = waux_mat(1:labm + 1, ima, ipbm - 1)
405 IF (ma > 0 .AND. mb - 1 > 0) wbmm(0:labm, 4) = waux_mat(1:labm + 1, ima, imbm + 1)
406 END IF
407 !*** Wbmp: lb-1, mb+1
408 IF (mb + 1 <= lb - 1) THEN
409 wbmp(0:labm, 1) = waux_mat(1:labm + 1, ipa, ipbm + 1)
410 IF (mb + 1 > 0) wbmp(0:labm, 2) = waux_mat(1:labm + 1, ipa, imbm - 1)
411 IF (ma > 0) wbmp(0:labm, 3) = waux_mat(1:labm + 1, ima, ipbm + 1)
412 IF (ma > 0 .AND. mb + 1 > 0) wbmp(0:labm, 4) = waux_mat(1:labm + 1, ima, imbm - 1)
413 END IF
414 DO j = 0, labmin
415 !*** x component
416 dwaux_mat(1, j + 1, ipa, ipb) = daa_p*wamp(j, 1) - daa_m*wamm(j, 1) &
417 - dab_p*wbmp(j, 1) + dab_m*wbmm(j, 1)
418 IF (mb > 0) THEN
419 dwaux_mat(1, j + 1, ipa, imb) = daa_p*wamp(j, 2) - daa_m*wamm(j, 2) &
420 - dab_p*wbmp(j, 2) + dab_m*wbmm(j, 2)
421 END IF
422 IF (ma > 0) THEN
423 dwaux_mat(1, j + 1, ima, ipb) = daa_p*wamp(j, 3) - daa_m*wamm(j, 3) &
424 - dab_p*wbmp(j, 3) + dab_m*wbmm(j, 3)
425 END IF
426 IF (ma > 0 .AND. mb > 0) THEN
427 dwaux_mat(1, j + 1, ima, imb) = daa_p*wamp(j, 4) - daa_m*wamm(j, 4) &
428 - dab_p*wbmp(j, 4) + dab_m*wbmm(j, 4)
429 END IF
430
431 !**** y component
432 dwaux_mat(2, j + 1, ipa, ipb) = daa_p*wamp(j, 3) + daa_m*wamm(j, 3) &
433 - dab_p*wbmp(j, 2) - dab_m*wbmm(j, 2)
434 IF (mb > 0) THEN
435 dwaux_mat(2, j + 1, ipa, imb) = daa_p*wamp(j, 4) + daa_m*wamm(j, 4) &
436 + dab_p*wbmp(j, 1) + dab_m*wbmm(j, 1)
437 END IF
438 IF (ma > 0) THEN
439 dwaux_mat(2, j + 1, ima, ipb) = -daa_p*wamp(j, 1) - daa_m*wamm(j, 1) &
440 - dab_p*wbmp(j, 4) - dab_m*wbmm(j, 4)
441 END IF
442 IF (ma > 0 .AND. mb > 0) THEN
443 dwaux_mat(2, j + 1, ima, imb) = -daa_p*wamp(j, 2) - daa_m*wamm(j, 2) &
444 + dab_p*wbmp(j, 3) + dab_m*wbmm(j, 3)
445 END IF
446 !**** z compnent
447 dwaux_mat(3, j + 1, ipa, ipb) = daa*wam(j, 1) - dab*wbm(j, 1)
448 IF (mb > 0) THEN
449 dwaux_mat(3, j + 1, ipa, imb) = daa*wam(j, 2) - dab*wbm(j, 2)
450 END IF
451 IF (ma > 0) THEN
452 dwaux_mat(3, j + 1, ima, ipb) = daa*wam(j, 3) - dab*wbm(j, 3)
453 END IF
454 IF (ma > 0 .AND. mb > 0) THEN
455 dwaux_mat(3, j + 1, ima, imb) = daa*wam(j, 4) - dab*wbm(j, 4)
456 END IF
457
458 END DO
459 END DO
460 END DO
461 END DO
462 END DO
463
464 DEALLOCATE (wam, wamm, wamp)
465 DEALLOCATE (wbm, wbmm, wbmp)
466 DEALLOCATE (da, da_p, da_m)
467
468 END SUBROUTINE get_dw_matrix
469
470! **************************************************************************************************
471!> \brief calculates [ab] SHG overlap integrals using precomputed angular-
472!> dependent part
473!> \param la set of l quantum number on a
474!> \param first_sgfa indexing
475!> \param nshella number of shells for a
476!> \param lb set of l quantum number on b
477!> \param first_sgfb indexing
478!> \param nshellb number of shells for b
479!> \param swork_cont contracted and normalized [s|s] integrals
480!> \param Waux_mat precomputed angular-dependent part
481!> \param sab contracted integral of spherical harmonic Gaussianslm
482! **************************************************************************************************
483 SUBROUTINE construct_int_shg_ab(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
484 swork_cont, Waux_mat, sab)
485
486 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
487 INTEGER, INTENT(IN) :: nshella
488 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
489 INTEGER, INTENT(IN) :: nshellb
490 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: swork_cont, waux_mat
491 REAL(kind=dp), DIMENSION(:, :), INTENT(INOUT) :: sab
492
493 INTEGER :: fnla, fnlb, fsgfa, fsgfb, ishella, j, &
494 jshellb, labmin, lai, lbj, lnla, lnlb, &
495 lsgfa, lsgfb, mai, mbj
496 REAL(kind=dp) :: prefac
497
498 DO jshellb = 1, nshellb
499 lbj = lb(jshellb)
500 fnlb = nsoset(lbj - 1) + 1
501 lnlb = nsoset(lbj)
502 fsgfb = first_sgfb(jshellb)
503 lsgfb = fsgfb + 2*lbj
504 DO ishella = 1, nshella
505 lai = la(ishella)
506 fnla = nsoset(lai - 1) + 1
507 lnla = nsoset(lai)
508 fsgfa = first_sgfa(ishella)
509 lsgfa = fsgfa + 2*lai
510 labmin = min(lai, lbj)
511 DO mbj = 0, 2*lbj
512 DO mai = 0, 2*lai
513 DO j = 0, labmin
514 prefac = swork_cont(lai + lbj - j + 1, ishella, jshellb)
515 sab(fsgfa + mai, fsgfb + mbj) = sab(fsgfa + mai, fsgfb + mbj) &
516 + prefac*waux_mat(j + 1, fnla + mai, fnlb + mbj)
517 END DO
518 END DO
519 END DO
520 END DO
521 END DO
522
523 END SUBROUTINE construct_int_shg_ab
524
525! **************************************************************************************************
526!> \brief calculates derivatives of [ab] SHG overlap integrals using precomputed
527!> angular-dependent part
528!> \param la set of l quantum number on a
529!> \param first_sgfa indexing
530!> \param nshella number of shells for a
531!> \param lb set of l quantum number on b
532!> \param first_sgfb indexing
533!> \param nshellb number of shells for b
534!> \param rab distance vector Ra-Rb
535!> \param swork_cont contracted and normalized [s|s] integrals
536!> \param Waux_mat precomputed angular-dependent part
537!> \param dWaux_mat ...
538!> \param dsab derivative of contracted integral of spherical harmonic Gaussians
539! **************************************************************************************************
540 SUBROUTINE construct_dev_shg_ab(la, first_sgfa, nshella, lb, first_sgfb, nshellb, rab, &
541 swork_cont, Waux_mat, dWaux_mat, dsab)
542
543 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
544 INTEGER, INTENT(IN) :: nshella
545 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
546 INTEGER, INTENT(IN) :: nshellb
547 REAL(kind=dp), INTENT(IN) :: rab(3)
548 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: swork_cont, waux_mat
549 REAL(kind=dp), DIMENSION(:, :, :, :), INTENT(IN) :: dwaux_mat
550 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: dsab
551
552 INTEGER :: fnla, fnlb, fsgfa, fsgfb, i, ishella, j, &
553 jshellb, labmin, lai, lbj, lnla, lnlb, &
554 lsgfa, lsgfb
555 REAL(kind=dp) :: dprefac, prefac, rabx2(3)
556
557 rabx2(:) = 2.0_dp*rab
558 DO jshellb = 1, nshellb
559 lbj = lb(jshellb)
560 fnlb = nsoset(lbj - 1) + 1
561 lnlb = nsoset(lbj)
562 fsgfb = first_sgfb(jshellb)
563 lsgfb = fsgfb + 2*lbj
564 DO ishella = 1, nshella
565 lai = la(ishella)
566 fnla = nsoset(lai - 1) + 1
567 lnla = nsoset(lai)
568 fsgfa = first_sgfa(ishella)
569 lsgfa = fsgfa + 2*lai
570 labmin = min(lai, lbj)
571 DO j = 0, labmin
572 prefac = swork_cont(lai + lbj - j + 1, ishella, jshellb)
573 dprefac = swork_cont(lai + lbj - j + 2, ishella, jshellb) !j+1
574 DO i = 1, 3
575 dsab(fsgfa:lsgfa, fsgfb:lsgfb, i) = dsab(fsgfa:lsgfa, fsgfb:lsgfb, i) &
576 + rabx2(i)*dprefac*waux_mat(j + 1, fnla:lnla, fnlb:lnlb) &
577 + prefac*dwaux_mat(i, j + 1, fnla:lnla, fnlb:lnlb)
578 END DO
579 END DO
580 END DO
581 END DO
582
583 END SUBROUTINE construct_dev_shg_ab
584
585! **************************************************************************************************
586!> \brief calculates [aba] SHG overlap integrals using precomputed angular-
587!> dependent part
588!> \param la set of l quantum number on a, orbital basis
589!> \param first_sgfa indexing
590!> \param nshella number of shells for a, orbital basis
591!> \param lb set of l quantum number on b. orbital basis
592!> \param first_sgfb indexing
593!> \param nshellb number of shells for b, orbital basis
594!> \param lca of l quantum number on a, aux basis
595!> \param first_sgfca indexing
596!> \param nshellca number of shells for a, aux basis
597!> \param cg_coeff Clebsch-Gordon coefficients
598!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
599!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
600!> \param swork_cont contracted and normalized [s|ra^n|s] integrals
601!> \param Waux_mat precomputed angular-dependent part
602!> \param saba contracted overlap [aba] of spherical harmonic Gaussians
603! **************************************************************************************************
604 SUBROUTINE construct_overlap_shg_aba(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
605 lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, &
606 ncg_none0, swork_cont, Waux_mat, saba)
607
608 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
609 INTEGER, INTENT(IN) :: nshella
610 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
611 INTEGER, INTENT(IN) :: nshellb
612 INTEGER, DIMENSION(:), INTENT(IN) :: lca, first_sgfca
613 INTEGER, INTENT(IN) :: nshellca
614 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
615 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
616 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
617 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
618 INTENT(IN) :: swork_cont
619 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
620 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: saba
621
622 INTEGER :: ia, il, ilist, ishella, isoa1, isoa2, isoaa, j, jb, jshellb, ka, kshella, laa, &
623 labmin, lai, lak, lbj, maa, mai, mak, mbj, nla, nlb, sgfa, sgfb, sgfca
624 REAL(kind=dp) :: prefac, stemp
625
626 DO kshella = 1, nshellca
627 lak = lca(kshella)
628 sgfca = first_sgfca(kshella)
629 ka = sgfca + lak
630 DO jshellb = 1, nshellb
631 lbj = lb(jshellb)
632 nlb = nsoset(lbj - 1) + lbj + 1
633 sgfb = first_sgfb(jshellb)
634 jb = sgfb + lbj
635 DO ishella = 1, nshella
636 lai = la(ishella)
637 sgfa = first_sgfa(ishella)
638 ia = sgfa + lai
639 DO mai = -lai, lai, 1
640 DO mak = -lak, lak, 1
641 isoa1 = indso_inv(lai, mai)
642 isoa2 = indso_inv(lak, mak)
643 DO mbj = -lbj, lbj, 1
644 DO ilist = 1, ncg_none0(isoa1, isoa2)
645 isoaa = cg_none0_list(isoa1, isoa2, ilist)
646 laa = indso(1, isoaa)
647 maa = indso(2, isoaa)
648 nla = nsoset(laa - 1) + laa + 1
649 labmin = min(laa, lbj)
650 il = int((lai + lak - laa)/2)
651 stemp = 0.0_dp
652 DO j = 0, labmin
653 prefac = swork_cont(laa + lbj - j + 1, il, ishella, jshellb, kshella)
654 stemp = stemp + prefac*waux_mat(j + 1, nla + maa, nlb + mbj)
655 END DO
656 saba(ia + mai, jb + mbj, ka + mak) = saba(ia + mai, jb + mbj, ka + mak) + cg_coeff(isoa1, isoa2, isoaa)*stemp
657 END DO
658 END DO
659 END DO
660 END DO
661 END DO
662 END DO
663 END DO
664
665 END SUBROUTINE construct_overlap_shg_aba
666
667! **************************************************************************************************
668!> \brief calculates derivatives of [aba] SHG overlap integrals using
669!> precomputed angular-dependent part
670!> \param la set of l quantum number on a, orbital basis
671!> \param first_sgfa indexing
672!> \param nshella number of shells for a, orbital basis
673!> \param lb set of l quantum number on b. orbital basis
674!> \param first_sgfb indexing
675!> \param nshellb number of shells for b, orbital basis
676!> \param lca of l quantum number on a, aux basis
677!> \param first_sgfca indexing
678!> \param nshellca number of shells for a, aux basis
679!> \param cg_coeff Clebsch-Gordon coefficients
680!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
681!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
682!> \param rab distance vector Ra-Rb
683!> \param swork_cont contracted and normalized [s|ra^n|s] integrals
684!> \param Waux_mat precomputed angular-dependent part
685!> \param dWaux_mat derivatives of precomputed angular-dependent part
686!> \param dsaba derivative of contracted overlap [aba] of spherical harmonic
687!> Gaussians
688! **************************************************************************************************
689 SUBROUTINE dev_overlap_shg_aba(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
690 lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, &
691 ncg_none0, rab, swork_cont, Waux_mat, dWaux_mat, dsaba)
692
693 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
694 INTEGER, INTENT(IN) :: nshella
695 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
696 INTEGER, INTENT(IN) :: nshellb
697 INTEGER, DIMENSION(:), INTENT(IN) :: lca, first_sgfca
698 INTEGER, INTENT(IN) :: nshellca
699 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
700 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
701 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
702 REAL(kind=dp), INTENT(IN) :: rab(3)
703 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
704 INTENT(IN) :: swork_cont
705 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
706 REAL(kind=dp), DIMENSION(:, :, :, :), INTENT(IN) :: dwaux_mat
707 REAL(kind=dp), DIMENSION(:, :, :, :), &
708 INTENT(INOUT) :: dsaba
709
710 INTEGER :: i, ia, il, ilist, ishella, isoa1, isoa2, isoaa, j, jb, jshellb, ka, kshella, laa, &
711 labmin, lai, lak, lbj, maa, mai, mak, mbj, nla, nlb, sgfa, sgfb, sgfca
712 REAL(kind=dp) :: dprefac, dtemp(3), prefac, rabx2(3)
713
714 rabx2(:) = 2.0_dp*rab
715
716 DO kshella = 1, nshellca
717 lak = lca(kshella)
718 sgfca = first_sgfca(kshella)
719 ka = sgfca + lak
720 DO jshellb = 1, nshellb
721 lbj = lb(jshellb)
722 nlb = nsoset(lbj - 1) + lbj + 1
723 sgfb = first_sgfb(jshellb)
724 jb = sgfb + lbj
725 DO ishella = 1, nshella
726 lai = la(ishella)
727 sgfa = first_sgfa(ishella)
728 ia = sgfa + lai
729 DO mai = -lai, lai, 1
730 DO mak = -lak, lak, 1
731 isoa1 = indso_inv(lai, mai)
732 isoa2 = indso_inv(lak, mak)
733 DO mbj = -lbj, lbj, 1
734 DO ilist = 1, ncg_none0(isoa1, isoa2)
735 isoaa = cg_none0_list(isoa1, isoa2, ilist)
736 laa = indso(1, isoaa)
737 maa = indso(2, isoaa)
738 nla = nsoset(laa - 1) + laa + 1
739 labmin = min(laa, lbj)
740 il = (lai + lak - laa)/2 ! lai+lak-laa always even
741 dtemp = 0.0_dp
742 DO j = 0, labmin
743 prefac = swork_cont(laa + lbj - j + 1, il, ishella, jshellb, kshella)
744 dprefac = swork_cont(laa + lbj - j + 2, il, ishella, jshellb, kshella)
745 DO i = 1, 3
746 dtemp(i) = dtemp(i) + rabx2(i)*dprefac*waux_mat(j + 1, nla + maa, nlb + mbj) &
747 + prefac*dwaux_mat(i, j + 1, nla + maa, nlb + mbj)
748 END DO
749 END DO
750 DO i = 1, 3
751 dsaba(ia + mai, jb + mbj, ka + mak, i) = dsaba(ia + mai, jb + mbj, ka + mak, i) &
752 + cg_coeff(isoa1, isoa2, isoaa)*dtemp(i)
753 END DO
754 END DO
755 END DO
756 END DO
757 END DO
758 END DO
759 END DO
760 END DO
761
762 END SUBROUTINE dev_overlap_shg_aba
763
764! **************************************************************************************************
765!> \brief calculates [abb] SHG overlap integrals using precomputed angular-
766!> dependent part
767!> \param la set of l quantum number on a, orbital basis
768!> \param first_sgfa indexing
769!> \param nshella number of shells for a, orbital basis
770!> \param lb set of l quantum number on b. orbital basis
771!> \param first_sgfb indexing
772!> \param nshellb number of shells for b, orbital basis
773!> \param lcb l quantum number on b, aux basis
774!> \param first_sgfcb indexing
775!> \param nshellcb number of shells for b, aux basis
776!> \param cg_coeff Clebsch-Gordon coefficients
777!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
778!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
779!> \param swork_cont contracted and normalized [s|rb^n|s] integrals
780!> \param Waux_mat precomputed angular-dependent part
781!> \param sabb contracted overlap [abb] of spherical harmonic Gaussians
782! **************************************************************************************************
783 SUBROUTINE construct_overlap_shg_abb(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
784 lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, &
785 ncg_none0, swork_cont, Waux_mat, sabb)
786
787 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
788 INTEGER, INTENT(IN) :: nshella
789 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
790 INTEGER, INTENT(IN) :: nshellb
791 INTEGER, DIMENSION(:), INTENT(IN) :: lcb, first_sgfcb
792 INTEGER, INTENT(IN) :: nshellcb
793 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
794 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
795 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
796 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
797 INTENT(IN) :: swork_cont
798 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
799 REAL(kind=dp), DIMENSION(:, :, :), INTENT(INOUT) :: sabb
800
801 INTEGER :: ia, il, ilist, ishella, isob1, isob2, isobb, j, jb, jshellb, kb, kshellb, labmin, &
802 lai, lbb, lbj, lbk, mai, mbb, mbj, mbk, nla, nlb, sgfa, sgfb, sgfcb
803 REAL(kind=dp) :: prefac, stemp, tsign
804
805 DO kshellb = 1, nshellcb
806 lbk = lcb(kshellb)
807 sgfcb = first_sgfcb(kshellb)
808 kb = sgfcb + lbk
809 DO jshellb = 1, nshellb
810 lbj = lb(jshellb)
811 sgfb = first_sgfb(jshellb)
812 jb = sgfb + lbj
813 DO ishella = 1, nshella
814 lai = la(ishella)
815 nla = nsoset(lai - 1) + lai + 1
816 sgfa = first_sgfa(ishella)
817 ia = sgfa + lai
818 DO mbj = -lbj, lbj, 1
819 DO mbk = -lbk, lbk, 1
820 isob1 = indso_inv(lbj, mbj)
821 isob2 = indso_inv(lbk, mbk)
822 DO mai = -lai, lai, 1
823 DO ilist = 1, ncg_none0(isob1, isob2)
824 isobb = cg_none0_list(isob1, isob2, ilist)
825 lbb = indso(1, isobb)
826 mbb = indso(2, isobb)
827 nlb = nsoset(lbb - 1) + lbb + 1
828 ! tsgin: because we take the transpose of auxmat (calculated for (la,lb))
829 tsign = 1.0_dp
830 IF (modulo(lbb - lai, 2) /= 0) tsign = -1.0_dp
831 labmin = min(lai, lbb)
832 il = int((lbj + lbk - lbb)/2)
833 stemp = 0.0_dp
834 DO j = 0, labmin
835 prefac = swork_cont(lai + lbb - j + 1, il, ishella, jshellb, kshellb)
836 stemp = stemp + prefac*waux_mat(j + 1, nlb + mbb, nla + mai)
837 END DO
838 sabb(ia + mai, jb + mbj, kb + mbk) = sabb(ia + mai, jb + mbj, kb + mbk) + tsign*cg_coeff(isob1, isob2, isobb)*stemp
839 END DO
840 END DO
841 END DO
842 END DO
843 END DO
844 END DO
845 END DO
846
847 END SUBROUTINE construct_overlap_shg_abb
848
849! **************************************************************************************************
850!> \brief calculates derivatives of [abb] SHG overlap integrals using
851!> precomputed angular-dependent part
852!> \param la set of l quantum number on a, orbital basis
853!> \param first_sgfa indexing
854!> \param nshella number of shells for a, orbital basis
855!> \param lb set of l quantum number on b. orbital basis
856!> \param first_sgfb indexing
857!> \param nshellb number of shells for b, orbital basis
858!> \param lcb l quantum number on b, aux basis
859!> \param first_sgfcb indexing
860!> \param nshellcb number of shells for b, aux basis
861!> \param cg_coeff Clebsch-Gordon coefficients
862!> \param cg_none0_list list of none-zero Clebsch-Gordon coefficients
863!> \param ncg_none0 number of non-zero Clebsch-Gordon coefficients
864!> \param rab distance vector Ra-Rb
865!> \param swork_cont contracted and normalized [s|rb^n|s] integrals
866!> \param Waux_mat precomputed angular-dependent part
867!> \param dWaux_mat derivatives of precomputed angular-dependent part
868!> \param dsabb derivative of contracted overlap [abb] of spherical harmonic
869!> Gaussians
870! **************************************************************************************************
871 SUBROUTINE dev_overlap_shg_abb(la, first_sgfa, nshella, lb, first_sgfb, nshellb, &
872 lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, &
873 ncg_none0, rab, swork_cont, Waux_mat, dWaux_mat, dsabb)
874
875 INTEGER, DIMENSION(:), INTENT(IN) :: la, first_sgfa
876 INTEGER, INTENT(IN) :: nshella
877 INTEGER, DIMENSION(:), INTENT(IN) :: lb, first_sgfb
878 INTEGER, INTENT(IN) :: nshellb
879 INTEGER, DIMENSION(:), INTENT(IN) :: lcb, first_sgfcb
880 INTEGER, INTENT(IN) :: nshellcb
881 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: cg_coeff
882 INTEGER, DIMENSION(:, :, :), INTENT(IN) :: cg_none0_list
883 INTEGER, DIMENSION(:, :), INTENT(IN) :: ncg_none0
884 REAL(kind=dp), INTENT(IN) :: rab(3)
885 REAL(kind=dp), DIMENSION(:, 0:, :, :, :), &
886 INTENT(IN) :: swork_cont
887 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: waux_mat
888 REAL(kind=dp), DIMENSION(:, :, :, :), INTENT(IN) :: dwaux_mat
889 REAL(kind=dp), DIMENSION(:, :, :, :), &
890 INTENT(INOUT) :: dsabb
891
892 INTEGER :: i, ia, il, ilist, ishella, isob1, isob2, isobb, j, jb, jshellb, kb, kshellb, &
893 labmin, lai, lbb, lbj, lbk, mai, mbb, mbj, mbk, nla, nlb, sgfa, sgfb, sgfcb
894 REAL(kind=dp) :: dprefac, dtemp(3), prefac, rabx2(3), &
895 tsign
896
897 rabx2(:) = 2.0_dp*rab
898
899 DO kshellb = 1, nshellcb
900 lbk = lcb(kshellb)
901 sgfcb = first_sgfcb(kshellb)
902 kb = sgfcb + lbk
903 DO jshellb = 1, nshellb
904 lbj = lb(jshellb)
905 sgfb = first_sgfb(jshellb)
906 jb = sgfb + lbj
907 DO ishella = 1, nshella
908 lai = la(ishella)
909 nla = nsoset(lai - 1) + lai + 1
910 sgfa = first_sgfa(ishella)
911 ia = sgfa + lai
912 DO mbj = -lbj, lbj, 1
913 DO mbk = -lbk, lbk, 1
914 isob1 = indso_inv(lbj, mbj)
915 isob2 = indso_inv(lbk, mbk)
916 DO mai = -lai, lai, 1
917 DO ilist = 1, ncg_none0(isob1, isob2)
918 isobb = cg_none0_list(isob1, isob2, ilist)
919 lbb = indso(1, isobb)
920 mbb = indso(2, isobb)
921 nlb = nsoset(lbb - 1) + lbb + 1
922 ! tsgin: because we take the transpose of auxmat (calculated for (la,lb))
923 tsign = 1.0_dp
924 IF (modulo(lbb - lai, 2) /= 0) tsign = -1.0_dp
925 labmin = min(lai, lbb)
926 il = (lbj + lbk - lbb)/2
927 dtemp = 0.0_dp
928 DO j = 0, labmin
929 prefac = swork_cont(lai + lbb - j + 1, il, ishella, jshellb, kshellb)
930 dprefac = swork_cont(lai + lbb - j + 2, il, ishella, jshellb, kshellb)
931 DO i = 1, 3
932 dtemp(i) = dtemp(i) + rabx2(i)*dprefac*waux_mat(j + 1, nlb + mbb, nla + mai) &
933 + prefac*dwaux_mat(i, j + 1, nlb + mbb, nla + mai)
934 END DO
935 END DO
936 DO i = 1, 3
937 dsabb(ia + mai, jb + mbj, kb + mbk, i) = dsabb(ia + mai, jb + mbj, kb + mbk, i) &
938 + tsign*cg_coeff(isob1, isob2, isobb)*dtemp(i)
939 END DO
940 END DO
941 END DO
942 END DO
943 END DO
944 END DO
945 END DO
946 END DO
947
948 END SUBROUTINE dev_overlap_shg_abb
949
950END MODULE construct_shg
951
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:117
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:692
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:542
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:786
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:485
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:874
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:607
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