26 #include "../base/base_uses.f90"
32 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'construct_shg'
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
57 INTEGER :: li, mi, prefac
58 REAL(kind=
dp) :: rc, rc_00, rlm, rmlm, rplm, rs, rs_00, &
70 rc = -0.5_dp*r(1)*rc_00
71 rs = -0.5_dp*r(2)*rc_00
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))
83 IF (
modulo(li, 2) /= 0)
THEN
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
99 rlm_c(mi + 1, -mi) = rlm
102 prefac = (li + mi)*(li - mi)
103 rplm = (real(2*li - 1,
dp)*r(3)*rlm - r2*rmlm)/real(prefac,
dp)
107 IF (
modulo(mi, 2) /= 0)
THEN
108 rlm_c(li, -mi) = -rlm
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
121 rlm_s(mi + 1, -mi) = -rlm
124 prefac = (li + mi)*(li - mi)
125 rplm = (real(2*li - 1,
dp)*r(3)*rlm - r2*rmlm)/real(prefac,
dp)
129 IF (
modulo(mi, 2) /= 0)
THEN
132 rlm_s(li, -mi) = -rlm
144 SUBROUTINE get_alm(lmax, A)
146 INTEGER,
INTENT(IN) :: lmax
147 REAL(kind=
dp),
DIMENSION(0:lmax, 0:lmax), &
151 REAL(kind=
dp) :: temp
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)
162 END SUBROUTINE get_alm
174 SUBROUTINE get_da_prefactors(lmax, dA_p, dA_m, dA)
176 INTEGER,
INTENT(IN) :: lmax
177 REAL(kind=
dp),
DIMENSION(0:lmax, 0:lmax), &
178 INTENT(INOUT) :: da_p, da_m, da
181 REAL(kind=
dp) :: bm, bm_m, bm_p
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)
197 END SUBROUTINE get_da_prefactors
212 INTEGER,
DIMENSION(:),
POINTER :: lamax
213 INTEGER,
INTENT(IN) :: lbmax, lmax
214 REAL(kind=
dp),
DIMENSION(0:lmax, -2*lmax:2*lmax), &
216 REAL(kind=
dp),
DIMENSION(:, :, :),
INTENT(INOUT) :: waux_mat
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
228 ALLOCATE (a(0:lmax, 0:lmax))
229 CALL get_alm(lmax, a)
238 IF (
modulo(lb, 2) /= 0) a_lbmb = -a_lbmb
241 alm_fac = a_lama*a_lbmb
245 prefac = alm_fac*real(2**(la + lb - j),
dp)*
dfac(2*j - 1)
251 IF (laj < abs(ma_m) .AND. laj < abs(ma_p)) cycle
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
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))
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))
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))
286 IF (ma > 0 .AND. mb > 0)
THEN
287 wmat(4) = wmat(4) + prefac/a_jk*(wa(3)*wb(3) + wa(4)*wb(4))
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)
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
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, &
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))
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)
342 IF (lb > 0) nlbm =
nsoset(lb - 2)
346 IF (la > 0) nlam =
nsoset(la - 2)
348 lamb = min(la - 1, lb)
349 labm = min(la, lb - 1)
354 ipb = nlb + lb + mb + 1
355 imb = nlb + lb - mb + 1
356 ipbm = nlbm + lb + mb
357 imbm = nlbm + lb - mb
362 ipa = nla + la + ma + 1
363 ima = nla + la - ma + 1
364 ipam = nlam + la + ma
365 imam = nlam + la - 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)
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)
381 IF (ma - 1 > 0 .AND. mb > 0) wamm(0:lamb, 4) = waux_mat(1:lamb + 1, imam + 1, imb)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
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)
447 dwaux_mat(3, j + 1, ipa, ipb) = daa*wam(j, 1) - dab*wbm(j, 1)
449 dwaux_mat(3, j + 1, ipa, imb) = daa*wam(j, 2) - dab*wbm(j, 2)
452 dwaux_mat(3, j + 1, ima, ipb) = daa*wam(j, 3) - dab*wbm(j, 3)
454 IF (ma > 0 .AND. mb > 0)
THEN
455 dwaux_mat(3, j + 1, ima, imb) = daa*wam(j, 4) - dab*wbm(j, 4)
464 DEALLOCATE (wam, wamm, wamp)
465 DEALLOCATE (wbm, wbmm, wbmp)
466 DEALLOCATE (da, da_p, da_m)
484 swork_cont, Waux_mat, sab)
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
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
498 DO jshellb = 1, nshellb
500 fnlb =
nsoset(lbj - 1) + 1
502 fsgfb = first_sgfb(jshellb)
503 lsgfb = fsgfb + 2*lbj
504 DO ishella = 1, nshella
506 fnla =
nsoset(lai - 1) + 1
508 fsgfa = first_sgfa(ishella)
509 lsgfa = fsgfa + 2*lai
510 labmin = min(lai, lbj)
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)
541 swork_cont, Waux_mat, dWaux_mat, dsab)
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
552 INTEGER :: fnla, fnlb, fsgfa, fsgfb, i, ishella, j, &
553 jshellb, labmin, lai, lbj, lnla, lnlb, &
555 REAL(kind=
dp) :: dprefac, prefac, rabx2(3)
557 rabx2(:) = 2.0_dp*rab
558 DO jshellb = 1, nshellb
560 fnlb =
nsoset(lbj - 1) + 1
562 fsgfb = first_sgfb(jshellb)
563 lsgfb = fsgfb + 2*lbj
564 DO ishella = 1, nshella
566 fnla =
nsoset(lai - 1) + 1
568 fsgfa = first_sgfa(ishella)
569 lsgfa = fsgfa + 2*lai
570 labmin = min(lai, lbj)
572 prefac = swork_cont(lai + lbj - j + 1, ishella, jshellb)
573 dprefac = swork_cont(lai + lbj - j + 2, ishella, jshellb)
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)
605 lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, &
606 ncg_none0, swork_cont, Waux_mat, saba)
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
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
626 DO kshella = 1, nshellca
628 sgfca = first_sgfca(kshella)
630 DO jshellb = 1, nshellb
632 nlb =
nsoset(lbj - 1) + lbj + 1
633 sgfb = first_sgfb(jshellb)
635 DO ishella = 1, nshella
637 sgfa = first_sgfa(ishella)
639 DO mai = -lai, lai, 1
640 DO mak = -lak, lak, 1
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)
653 prefac = swork_cont(laa + lbj - j + 1, il, ishella, jshellb, kshella)
654 stemp = stemp + prefac*waux_mat(j + 1, nla + maa, nlb + mbj)
656 saba(ia + mai, jb + mbj, ka + mak) = saba(ia + mai, jb + mbj, ka + mak) + cg_coeff(isoa1, isoa2, isoaa)*stemp
690 lca, first_sgfca, nshellca, cg_coeff, cg_none0_list, &
691 ncg_none0, rab, swork_cont, Waux_mat, dWaux_mat, dsaba)
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
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)
714 rabx2(:) = 2.0_dp*rab
716 DO kshella = 1, nshellca
718 sgfca = first_sgfca(kshella)
720 DO jshellb = 1, nshellb
722 nlb =
nsoset(lbj - 1) + lbj + 1
723 sgfb = first_sgfb(jshellb)
725 DO ishella = 1, nshella
727 sgfa = first_sgfa(ishella)
729 DO mai = -lai, lai, 1
730 DO mak = -lak, lak, 1
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
743 prefac = swork_cont(laa + lbj - j + 1, il, ishella, jshellb, kshella)
744 dprefac = swork_cont(laa + lbj - j + 2, il, ishella, jshellb, kshella)
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)
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)
784 lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, &
785 ncg_none0, swork_cont, Waux_mat, sabb)
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
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
805 DO kshellb = 1, nshellcb
807 sgfcb = first_sgfcb(kshellb)
809 DO jshellb = 1, nshellb
811 sgfb = first_sgfb(jshellb)
813 DO ishella = 1, nshella
815 nla =
nsoset(lai - 1) + lai + 1
816 sgfa = first_sgfa(ishella)
818 DO mbj = -lbj, lbj, 1
819 DO mbk = -lbk, lbk, 1
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
830 IF (
modulo(lbb - lai, 2) /= 0) tsign = -1.0_dp
831 labmin = min(lai, lbb)
832 il = int((lbj + lbk - lbb)/2)
835 prefac = swork_cont(lai + lbb - j + 1, il, ishella, jshellb, kshellb)
836 stemp = stemp + prefac*waux_mat(j + 1, nlb + mbb, nla + mai)
838 sabb(ia + mai, jb + mbj, kb + mbk) = sabb(ia + mai, jb + mbj, kb + mbk) + tsign*cg_coeff(isob1, isob2, isobb)*stemp
872 lcb, first_sgfcb, nshellcb, cg_coeff, cg_none0_list, &
873 ncg_none0, rab, swork_cont, Waux_mat, dWaux_mat, dsabb)
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
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), &
897 rabx2(:) = 2.0_dp*rab
899 DO kshellb = 1, nshellcb
901 sgfcb = first_sgfcb(kshellb)
903 DO jshellb = 1, nshellb
905 sgfb = first_sgfb(jshellb)
907 DO ishella = 1, nshella
909 nla =
nsoset(lai - 1) + lai + 1
910 sgfa = first_sgfa(ishella)
912 DO mbj = -lbj, lbj, 1
913 DO mbk = -lbk, lbk, 1
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
924 IF (
modulo(lbb - lai, 2) /= 0) tsign = -1.0_dp
925 labmin = min(lai, lbb)
926 il = (lbj + lbk - lbb)/2
929 prefac = swork_cont(lai + lbb - j + 1, il, ishella, jshellb, kshellb)
930 dprefac = swork_cont(lai + lbb - j + 2, il, ishella, jshellb, kshellb)
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)
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)
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
Calculation of the integrals over solid harmonic Gaussian(SHG) functions. Routines for (a|O(r12)|b) a...
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
subroutine, public get_w_matrix(lamax, lbmax, lmax, Rc, Rs, Waux_mat)
calculates the angular dependent-part of the SHG integrals, transformation matrix W,...
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
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
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
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
subroutine, public get_dw_matrix(lamax, lbmax, Waux_mat, dWaux_mat)
calculates derivatives of transformation matrix W,
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
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
Defines the basic variable types.
integer, parameter, public dp
Definition of mathematical constants and functions.
real(kind=dp), dimension(-1:2 *maxfac+1), parameter, public dfac
real(kind=dp), dimension(0:maxfac), parameter, public fac
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public nsoset
integer, dimension(:, :), allocatable, public indso
integer, dimension(:, :), allocatable, public indso_inv