(git:374b731)
Loading...
Searching...
No Matches
debug_os_integrals.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 Debugs Obara-Saika integral matrices
10!> \par History
11!> created [07.2014]
12!> \authors Dorothea Golze
13! **************************************************************************************************
14MODULE debug_os_integrals
15
16 USE ai_overlap, ONLY: overlap
17 USE ai_overlap3, ONLY: overlap3
18 USE ai_overlap3_debug, ONLY: init_os_overlap3,&
19 os_overlap3
20 USE ai_overlap_aabb, ONLY: overlap_aabb
21 USE ai_overlap_debug, ONLY: init_os_overlap2,&
22 os_overlap2
23 USE kinds, ONLY: dp
24 USE orbital_pointers, ONLY: coset,&
25 indco,&
26 ncoset
27#include "./base/base_uses.f90"
28
29 IMPLICIT NONE
30
31 PRIVATE
32
33! **************************************************************************************************
34
35 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'debug_os_integrals'
36
37 PUBLIC :: overlap_ab_test, overlap_abc_test, overlap_aabb_test
38
39! **************************************************************************************************
40
41CONTAINS
42
43! ***************************************************************************************************
44!> \brief recursive test routines for integral (a,b) for only two exponents
45! **************************************************************************************************
46 SUBROUTINE overlap_ab_test_simple()
47
48 INTEGER :: ia1, iax, iay, iaz, ib1, ibx, iby, ibz, &
49 la_max, la_min, lb_max, lb_min, lds, &
50 ma, maxl, mb
51 INTEGER, DIMENSION(3) :: na, nb
52 REAL(KIND=dp) :: dab, dmax, res1, xa, xb
53 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: sab
54 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: swork
55 REAL(KIND=dp), DIMENSION(1) :: rpgfa, rpgfb, xa_work, xb_work
56 REAL(KIND=dp), DIMENSION(3) :: a, b, rab
57
58 xa = 0.783300000000_dp ! exponents
59 xb = 1.239648746700_dp
60
61 a = (/0.329309000000_dp, 0.28408240000_dp, 0.28408240000_dp/) !* bohr !positions
62 b = (/0.983983000000_dp, 0.00453720000_dp, 0.00432740000_dp/) !* bohr
63
64 la_min = 0
65 lb_min = 0
66
67 la_max = 3
68 lb_max = 4
69
70 !---------------------------------------
71 ALLOCATE (sab(ncoset(la_max), ncoset(lb_max)))
72
73 maxl = max(la_max, lb_max)
74 lds = ncoset(maxl)
75 ALLOCATE (swork(lds, lds, 1))
76 sab = 0._dp
77 rab(:) = b(:) - a(:)
78 dab = sqrt(dot_product(rab, rab))
79 xa_work(1) = xa
80 xb_work(1) = xb
81 rpgfa = 20._dp
82 rpgfb = 20._dp
83 CALL overlap(la_max_set=la_max, la_min_set=la_min, npgfa=1, rpgfa=rpgfa, zeta=xa_work, &
84 lb_max_set=lb_max, lb_min_set=lb_min, npgfb=1, rpgfb=rpgfb, zetb=xb_work, &
85 rab=rab, dab=dab, sab=sab, da_max_set=0, return_derivatives=.false., s=swork, lds=lds)
86 !---------------------------------------
87
88 CALL init_os_overlap2(xa, xb, a, b)
89
90 dmax = 0._dp
91 DO ma = la_min, la_max
92 DO mb = lb_min, lb_max
93 DO iax = 0, ma
94 DO iay = 0, ma - iax
95 iaz = ma - iax - iay
96 na(1) = iax; na(2) = iay; na(3) = iaz
97 ia1 = coset(iax, iay, iaz)
98 DO ibx = 0, mb
99 DO iby = 0, mb - ibx
100 ibz = mb - ibx - iby
101 nb(1) = ibx; nb(2) = iby; nb(3) = ibz
102 ib1 = coset(ibx, iby, ibz)
103 res1 = os_overlap2(na, nb)
104 ! write(*,*) "la, lb,na, nb, res1", ma, mb, na, nb, res1
105 ! write(*,*) "sab ia1, ib1", ia1, ib1, sab(ia1,ib1)
106 dmax = max(dmax, abs(res1 - sab(ia1, ib1)))
107 END DO
108 END DO
109 END DO
110 END DO
111 END DO
112 END DO
113
114 DEALLOCATE (sab, swork)
115
116 END SUBROUTINE overlap_ab_test_simple
117
118! ***************************************************************************************************
119!> \brief recursive test routines for integral (a,b)
120!> \param la_max ...
121!> \param la_min ...
122!> \param npgfa ...
123!> \param zeta ...
124!> \param lb_max ...
125!> \param lb_min ...
126!> \param npgfb ...
127!> \param zetb ...
128!> \param ra ...
129!> \param rb ...
130!> \param sab ...
131!> \param dmax ...
132! **************************************************************************************************
133 SUBROUTINE overlap_ab_test(la_max, la_min, npgfa, zeta, lb_max, lb_min, npgfb, zetb, &
134 ra, rb, sab, dmax)
135
136 INTEGER, INTENT(IN) :: la_max, la_min, npgfa
137 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zeta
138 INTEGER, INTENT(IN) :: lb_max, lb_min, npgfb
139 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zetb
140 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: ra, rb
141 REAL(kind=dp), DIMENSION(:, :), INTENT(IN) :: sab
142 REAL(kind=dp), INTENT(INOUT) :: dmax
143
144 INTEGER :: coa, cob, ia1, iax, iay, iaz, ib1, ibx, &
145 iby, ibz, ipgf, jpgf, ma, mb
146 INTEGER, DIMENSION(3) :: na, nb
147 REAL(kind=dp) :: res1, res2, xa, xb
148 REAL(kind=dp), DIMENSION(3) :: a, b
149
150 coa = 0
151 DO ipgf = 1, npgfa
152 cob = 0
153 DO jpgf = 1, npgfb
154 xa = zeta(ipgf) !exponents
155 xb = zetb(jpgf)
156 a = ra !positions
157 b = rb
158 CALL init_os_overlap2(xa, xb, a, b)
159 DO ma = la_min, la_max
160 DO mb = lb_min, lb_max
161 DO iax = 0, ma
162 DO iay = 0, ma - iax
163 iaz = ma - iax - iay
164 na(1) = iax; na(2) = iay; na(3) = iaz
165 ia1 = coset(iax, iay, iaz)
166 DO ibx = 0, mb
167 DO iby = 0, mb - ibx
168 ibz = mb - ibx - iby
169 nb(1) = ibx; nb(2) = iby; nb(3) = ibz
170 ib1 = coset(ibx, iby, ibz)
171 res1 = os_overlap2(na, nb)
172 res2 = sab(coa + ia1, cob + ib1)
173 dmax = max(dmax, abs(res1 - res2))
174 END DO
175 END DO
176 END DO
177 END DO
178 END DO
179 END DO
180 cob = cob + ncoset(lb_max)
181 END DO
182 coa = coa + ncoset(la_max)
183 END DO
184 !WRITE(*,*) "dmax overlap_ab_test", dmax
185
186 END SUBROUTINE overlap_ab_test
187
188! ***************************************************************************************************
189!> \brief recursive test routines for integral (a,b,c) for only three exponents
190! **************************************************************************************************
191 SUBROUTINE overlap_abc_test_simple()
192
193 INTEGER :: ia1, iax, iay, iaz, ib1, ibx, iby, ibz, &
194 ic1, icx, icy, icz, la_max, la_min, &
195 lb_max, lb_min, lc_max, lc_min, ma, &
196 mb, mc
197 INTEGER, DIMENSION(3) :: na, nb, nc
198 REAL(kind=dp) :: dab, dac, dbc, dmax, res1, xa, xb, xc
199 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :) :: sabc
200 REAL(kind=dp), DIMENSION(1) :: rpgfa, rpgfb, rpgfc, xa_work, xb_work, &
201 xc_work
202 REAL(kind=dp), DIMENSION(3) :: a, b, c, rab, rac, rbc
203
204 xa = 0.783300000000_dp ! exponents
205 xb = 1.239648746700_dp
206 xc = 0.548370000000_dp
207
208 a = (/0.329309000000_dp, 0.28408240000_dp, 0.28408240000_dp/) !* bohr !positions
209 b = (/0.983983000000_dp, 0.00453720000_dp, 0.00432740000_dp/) !* bohr
210 c = (/0.032380000000_dp, 1.23470000000_dp, 0.11137400000_dp/) !* bohr
211
212 la_min = 0
213 lb_min = 0
214 lc_min = 0
215
216 la_max = 0
217 lb_max = 0
218 lc_max = 1
219
220 !---------------------------------------
221 rab(:) = b(:) - a(:)
222 dab = sqrt(dot_product(rab, rab))
223 rac(:) = c(:) - a(:)
224 dac = sqrt(dot_product(rac, rac))
225 rbc(:) = c(:) - b(:)
226 dbc = sqrt(dot_product(rbc, rbc))
227 ALLOCATE (sabc(ncoset(la_max), ncoset(lb_max), ncoset(lc_max)))
228 xa_work(1) = xa
229 xb_work(1) = xb
230 xc_work(1) = xc
231 rpgfa = 20._dp
232 rpgfb = 20._dp
233 rpgfc = 20._dp
234 sabc = 0._dp
235
236 CALL overlap3(la_max_set=la_max, npgfa=1, zeta=xa_work, rpgfa=rpgfa, la_min_set=la_min, &
237 lb_max_set=lb_max, npgfb=1, zetb=xb_work, rpgfb=rpgfb, lb_min_set=lb_min, &
238 lc_max_set=lc_max, npgfc=1, zetc=xc_work, rpgfc=rpgfc, lc_min_set=lc_min, &
239 rab=rab, dab=dab, rac=rac, dac=dac, rbc=rbc, dbc=dbc, sabc=sabc)
240
241 !---------------------------------------
242
243 CALL init_os_overlap3(xa, xb, xc, a, b, c)
244
245 dmax = 0._dp
246 DO ma = la_min, la_max
247 DO mc = lc_min, lc_max
248 DO mb = lb_min, lb_max
249 DO iax = 0, ma
250 DO iay = 0, ma - iax
251 iaz = ma - iax - iay
252 na(1) = iax; na(2) = iay; na(3) = iaz
253 ia1 = coset(iax, iay, iaz)
254 DO icx = 0, mc
255 DO icy = 0, mc - icx
256 icz = mc - icx - icy
257 nc(1) = icx; nc(2) = icy; nc(3) = icz
258 ic1 = coset(icx, icy, icz)
259 DO ibx = 0, mb
260 DO iby = 0, mb - ibx
261 ibz = mb - ibx - iby
262 nb(1) = ibx; nb(2) = iby; nb(3) = ibz
263 ib1 = coset(ibx, iby, ibz)
264 res1 = os_overlap3(na, nc, nb)
265 !write(*,*) "la, lc, lb,na,nc, nb, res1", ma, mc, mb, na, nc, nb, res1
266 !write(*,*) "sabc ia1, ib1, ic1", ia1, ib1, ic1, sabc(ia1,ib1,ic1)
267 dmax = max(dmax, abs(res1 - sabc(ia1, ib1, ic1)))
268 END DO
269 END DO
270 END DO
271 END DO
272 END DO
273 END DO
274 END DO
275 END DO
276 END DO
277
278 DEALLOCATE (sabc)
279
280 END SUBROUTINE overlap_abc_test_simple
281
282! ***************************************************************************************************
283!> \brief recursive test routines for integral (a,b,c)
284!> \param la_max ...
285!> \param npgfa ...
286!> \param zeta ...
287!> \param la_min ...
288!> \param lb_max ...
289!> \param npgfb ...
290!> \param zetb ...
291!> \param lb_min ...
292!> \param lc_max ...
293!> \param npgfc ...
294!> \param zetc ...
295!> \param lc_min ...
296!> \param ra ...
297!> \param rb ...
298!> \param rc ...
299!> \param sabc ...
300!> \param dmax ...
301! **************************************************************************************************
302 SUBROUTINE overlap_abc_test(la_max, npgfa, zeta, la_min, &
303 lb_max, npgfb, zetb, lb_min, &
304 lc_max, npgfc, zetc, lc_min, &
305 ra, rb, rc, sabc, dmax)
306
307 INTEGER, INTENT(IN) :: la_max, npgfa
308 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zeta
309 INTEGER, INTENT(IN) :: la_min, lb_max, npgfb
310 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zetb
311 INTEGER, INTENT(IN) :: lb_min, lc_max, npgfc
312 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zetc
313 INTEGER, INTENT(IN) :: lc_min
314 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: ra, rb, rc
315 REAL(kind=dp), DIMENSION(:, :, :), INTENT(IN) :: sabc
316 REAL(kind=dp), INTENT(INOUT) :: dmax
317
318 INTEGER :: coa, cob, coc, ia1, iax, iay, iaz, ib1, &
319 ibx, iby, ibz, ic1, icx, icy, icz, &
320 ipgf, jpgf, kpgf, ma, mb, mc
321 INTEGER, DIMENSION(3) :: na, nb, nc
322 REAL(kind=dp) :: res1, res2, xa, xb, xc
323 REAL(kind=dp), DIMENSION(3) :: a, b, c
324
325 coa = 0
326 DO ipgf = 1, npgfa
327 cob = 0
328 DO jpgf = 1, npgfb
329 coc = 0
330 DO kpgf = 1, npgfc
331
332 xa = zeta(ipgf) ! exponents
333 xb = zetb(jpgf)
334 xc = zetc(kpgf)
335
336 a = ra !positions
337 b = rb
338 c = rc
339
340 CALL init_os_overlap3(xa, xb, xc, a, b, c)
341
342 DO ma = la_min, la_max
343 DO mc = lc_min, lc_max
344 DO mb = lb_min, lb_max
345 DO iax = 0, ma
346 DO iay = 0, ma - iax
347 iaz = ma - iax - iay
348 na(1) = iax; na(2) = iay; na(3) = iaz
349 ia1 = coset(iax, iay, iaz)
350 DO icx = 0, mc
351 DO icy = 0, mc - icx
352 icz = mc - icx - icy
353 nc(1) = icx; nc(2) = icy; nc(3) = icz
354 ic1 = coset(icx, icy, icz)
355 DO ibx = 0, mb
356 DO iby = 0, mb - ibx
357 ibz = mb - ibx - iby
358 nb(1) = ibx; nb(2) = iby; nb(3) = ibz
359 ib1 = coset(ibx, iby, ibz)
360 res1 = os_overlap3(na, nc, nb)
361 res2 = sabc(coa + ia1, cob + ib1, coc + ic1)
362 dmax = max(dmax, abs(res1 - res2))
363 !IF(dmax > 1.E-10) WRITE(*,*) "dmax in loop", dmax
364 END DO
365 END DO
366 END DO
367 END DO
368 END DO
369 END DO
370 END DO
371 END DO
372 END DO
373 coc = coc + ncoset(lc_max)
374 END DO
375 cob = cob + ncoset(lb_max)
376 END DO
377 coa = coa + ncoset(la_max)
378 END DO
379 !WRITE(*,*) "dmax abc", dmax
380
381 END SUBROUTINE overlap_abc_test
382
383! ***************************************************************************************************
384!> \brief recursive test routines for integral (aa,bb) for only four exponents
385! **************************************************************************************************
386 SUBROUTINE overlap_aabb_test_simple()
387
388 INTEGER :: i, iax, iay, iaz, ibx, iby, ibz, j, k, l, la_max, la_max1, la_max2, la_min, &
389 la_min1, la_min2, lb_max, lb_max1, lb_max2, lb_min, lb_min1, lb_min2, lds, ma, maxl, mb
390 INTEGER, DIMENSION(3) :: na, naa, nb, nbb
391 REAL(kind=dp) :: dab, dmax, res1, xa, xa1, xa2, xb, xb1, &
392 xb2
393 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: swork
394 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :, :) :: saabb
395 REAL(kind=dp), DIMENSION(1) :: rpgfa1, rpgfa2, rpgfb1, rpgfb2, &
396 xa_work1, xa_work2, xb_work1, xb_work2
397 REAL(kind=dp), DIMENSION(3) :: a, b, rab
398
399 xa1 = 0.783300000000_dp ! exponents
400 xb1 = 1.239648746700_dp
401 xa2 = 0.3400000000_dp ! exponents
402 xb2 = 2.76_dp
403
404 a = (/0.329309000000_dp, 0.28408240000_dp, 0.28408240000_dp/) !* bohr !positions
405 b = (/0.983983000000_dp, 0.00453720000_dp, 0.00432740000_dp/) !* bohr
406
407 la_min1 = 0
408 la_min2 = 0
409 lb_min1 = 3
410 lb_min2 = 1
411
412 la_max1 = 1
413 la_max2 = 2
414 lb_max1 = 3
415 lb_max2 = 4
416
417 !---------------------------------------
418 ALLOCATE (saabb(ncoset(la_max1), ncoset(la_max2), ncoset(lb_max1), ncoset(lb_max2)))
419
420 maxl = max(la_max1 + la_max2, lb_max1 + lb_max2)
421 lds = ncoset(maxl)
422 ALLOCATE (swork(lds, lds))
423 saabb = 0._dp
424 rab(:) = b(:) - a(:)
425 dab = sqrt(dot_product(rab, rab))
426 xa_work1(1) = xa1
427 xa_work2(1) = xa2
428 xb_work1(1) = xb1
429 xb_work2(1) = xb2
430 rpgfa1 = 20._dp
431 rpgfa2 = 20._dp
432 rpgfb1 = 20._dp
433 rpgfb2 = 20._dp
434 CALL overlap_aabb(la_max_set1=la_max1, la_min_set1=la_min1, npgfa1=1, rpgfa1=rpgfa1, zeta1=xa_work1, &
435 la_max_set2=la_max2, la_min_set2=la_min2, npgfa2=1, rpgfa2=rpgfa2, zeta2=xa_work2, &
436 lb_max_set1=lb_max1, lb_min_set1=lb_min1, npgfb1=1, rpgfb1=rpgfb1, zetb1=xb_work1, &
437 lb_max_set2=lb_max2, lb_min_set2=lb_min2, npgfb2=1, rpgfb2=rpgfb2, zetb2=xb_work2, &
438 asets_equal=.false., bsets_equal=.false., rab=rab, dab=dab, saabb=saabb, s=swork, lds=lds)
439 !---------------------------------------
440
441 xa = xa1 + xa2
442 xb = xb1 + xb2
443 la_min = la_min1 + la_min2
444 la_max = la_max1 + la_max2
445 lb_min = lb_min1 + lb_min2
446 lb_max = lb_max1 + lb_max2
447
448 CALL init_os_overlap2(xa, xb, a, b)
449
450 dmax = 0._dp
451 DO ma = la_min, la_max
452 DO mb = lb_min, lb_max
453 DO iax = 0, ma
454 DO iay = 0, ma - iax
455 iaz = ma - iax - iay
456 na(1) = iax; na(2) = iay; na(3) = iaz
457 DO ibx = 0, mb
458 DO iby = 0, mb - ibx
459 ibz = mb - ibx - iby
460 nb(1) = ibx; nb(2) = iby; nb(3) = ibz
461 res1 = os_overlap2(na, nb)
462 DO i = ncoset(la_min1 - 1) + 1, ncoset(la_max1)
463 DO j = ncoset(la_min2 - 1) + 1, ncoset(la_max2)
464 naa = indco(1:3, i) + indco(1:3, j)
465 DO k = ncoset(lb_min1 - 1) + 1, ncoset(lb_max1)
466 DO l = ncoset(lb_min2 - 1) + 1, ncoset(lb_max2)
467 nbb = indco(1:3, k) + indco(1:3, l)
468 IF (all(na == naa) .AND. all(nb == nbb)) THEN
469 dmax = max(dmax, abs(res1 - saabb(i, j, k, l)))
470 END IF
471 END DO
472 END DO
473 END DO
474 END DO
475 END DO
476 END DO
477 END DO
478 END DO
479 END DO
480 END DO
481
482 DEALLOCATE (saabb, swork)
483
484 END SUBROUTINE overlap_aabb_test_simple
485
486! ***************************************************************************************************
487!> \brief recursive test routines for integral (aa,bb)
488!> \param la_max1 ...
489!> \param la_min1 ...
490!> \param npgfa1 ...
491!> \param zeta1 ...
492!> \param la_max2 ...
493!> \param la_min2 ...
494!> \param npgfa2 ...
495!> \param zeta2 ...
496!> \param lb_max1 ...
497!> \param lb_min1 ...
498!> \param npgfb1 ...
499!> \param zetb1 ...
500!> \param lb_max2 ...
501!> \param lb_min2 ...
502!> \param npgfb2 ...
503!> \param zetb2 ...
504!> \param ra ...
505!> \param rb ...
506!> \param saabb ...
507!> \param dmax ...
508! **************************************************************************************************
509 SUBROUTINE overlap_aabb_test(la_max1, la_min1, npgfa1, zeta1, &
510 la_max2, la_min2, npgfa2, zeta2, &
511 lb_max1, lb_min1, npgfb1, zetb1, &
512 lb_max2, lb_min2, npgfb2, zetb2, &
513 ra, rb, saabb, dmax)
514
515 INTEGER, INTENT(IN) :: la_max1, la_min1, npgfa1
516 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zeta1
517 INTEGER, INTENT(IN) :: la_max2, la_min2, npgfa2
518 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zeta2
519 INTEGER, INTENT(IN) :: lb_max1, lb_min1, npgfb1
520 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zetb1
521 INTEGER, INTENT(IN) :: lb_max2, lb_min2, npgfb2
522 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zetb2
523 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: ra, rb
524 REAL(kind=dp), DIMENSION(:, :, :, :), INTENT(IN) :: saabb
525 REAL(kind=dp), INTENT(INOUT) :: dmax
526
527 INTEGER :: coa1, coa2, cob1, cob2, i, iax, iay, &
528 iaz, ibx, iby, ibz, ipgf, j, jpgf, k, &
529 kpgf, l, la_max, la_min, lb_max, &
530 lb_min, lpgf, ma, mb
531 INTEGER, DIMENSION(3) :: na, naa, nb, nbb
532 REAL(kind=dp) :: res1, xa, xb
533 REAL(kind=dp), DIMENSION(3) :: a, b
534
535 coa1 = 0
536 DO ipgf = 1, npgfa1
537 coa2 = 0
538 DO jpgf = 1, npgfa2
539 cob1 = 0
540 DO kpgf = 1, npgfb1
541 cob2 = 0
542 DO lpgf = 1, npgfb2
543
544 xa = zeta1(ipgf) + zeta2(jpgf) ! exponents
545 xb = zetb1(kpgf) + zetb2(lpgf) ! exponents
546 la_max = la_max1 + la_max2
547 lb_max = lb_max1 + lb_max2
548 la_min = la_min1 + la_min2
549 lb_min = lb_min1 + lb_min2
550
551 a = ra !positions
552 b = rb
553
554 CALL init_os_overlap2(xa, xb, a, b)
555
556 DO ma = la_min, la_max
557 DO mb = lb_min, lb_max
558 DO iax = 0, ma
559 DO iay = 0, ma - iax
560 iaz = ma - iax - iay
561 na(1) = iax; na(2) = iay; na(3) = iaz
562 DO ibx = 0, mb
563 DO iby = 0, mb - ibx
564 ibz = mb - ibx - iby
565 nb(1) = ibx; nb(2) = iby; nb(3) = ibz
566 res1 = os_overlap2(na, nb)
567 DO i = ncoset(la_min1 - 1) + 1, ncoset(la_max1)
568 DO j = ncoset(la_min2 - 1) + 1, ncoset(la_max2)
569 naa = indco(1:3, i) + indco(1:3, j)
570 DO k = ncoset(lb_min1 - 1) + 1, ncoset(lb_max1)
571 DO l = ncoset(lb_min2 - 1) + 1, ncoset(lb_max2)
572 nbb = indco(1:3, k) + indco(1:3, l)
573 IF (all(na == naa) .AND. all(nb == nbb)) THEN
574 dmax = max(dmax, abs(res1 - saabb(coa1 + i, coa2 + j, cob1 + k, cob2 + l)))
575 END IF
576 END DO
577 END DO
578 END DO
579 END DO
580 END DO
581 END DO
582 END DO
583 END DO
584 END DO
585 END DO
586 cob2 = cob2 + ncoset(lb_max2)
587 END DO
588 cob1 = cob1 + ncoset(lb_max1)
589 END DO
590 coa2 = coa2 + ncoset(la_max2)
591 END DO
592 coa1 = coa1 + ncoset(la_max1)
593 END DO
594
595 !WRITE(*,*) "dmax aabb", dmax
596
597 END SUBROUTINE overlap_aabb_test
598
599END MODULE debug_os_integrals
ai_overlap3_debug
Three-center integrals over Cartesian Gaussian-type functions.
Definition ai_overlap3_debug.F:16
ai_overlap3_debug::a
real(dp), dimension(3) a
Definition ai_overlap3_debug.F:29
ai_overlap3_debug::init_os_overlap3
subroutine, public init_os_overlap3(ya, yb, yc, ra, rb, rc)
Calculation of three-center integrals over Cartesian Gaussian-type functions.
Definition ai_overlap3_debug.F:49
ai_overlap3_debug::c
real(dp), dimension(3) c
Definition ai_overlap3_debug.F:29
ai_overlap3_debug::b
real(dp), dimension(3) b
Definition ai_overlap3_debug.F:29
ai_overlap3_debug::xb
real(dp) xb
Definition ai_overlap3_debug.F:28
ai_overlap3_debug::os_overlap3
recursive real(dp) function, public os_overlap3(an, cn, bn)
...
Definition ai_overlap3_debug.F:84
ai_overlap3_debug::xa
real(dp) xa
Definition ai_overlap3_debug.F:28
ai_overlap3_debug::zeta
real(dp) zeta
Definition ai_overlap3_debug.F:31
ai_overlap3
Definition ai_overlap3.F:24
ai_overlap3::overlap3
subroutine, public overlap3(la_max_set, npgfa, zeta, rpgfa, la_min_set, lb_max_set, npgfb, zetb, rpgfb, lb_min_set, lc_max_set, npgfc, zetc, rpgfc, lc_min_set, rab, dab, rac, dac, rbc, dbc, sabc, sdabc, sabdc, int_abc_ext)
Calculation of three-center overlap integrals [a|b|c] over primitive Cartesian Gaussian functions.
Definition ai_overlap3.F:107
ai_overlap_aabb
Calculation of the overlap integrals over Cartesian Gaussian-type functions.
Definition ai_overlap_aabb.F:27
ai_overlap_debug
Two-center overlap integrals over Cartesian Gaussian-type functions.
Definition ai_overlap_debug.F:16
ai_overlap_debug::init_os_overlap2
subroutine, public init_os_overlap2(ya, yb, ra, rb)
Calculation of overlap integrals over Cartesian Gaussian-type functions.
Definition ai_overlap_debug.F:47
ai_overlap_debug::os_overlap2
recursive real(dp) function, public os_overlap2(an, bn)
...
Definition ai_overlap_debug.F:73
ai_overlap
Calculation of the overlap integrals over Cartesian Gaussian-type functions.
Definition ai_overlap.F:18
ai_overlap::overlap
subroutine, public overlap(la_max_set, la_min_set, npgfa, rpgfa, zeta, lb_max_set, lb_min_set, npgfb, rpgfb, zetb, rab, dab, sab, da_max_set, return_derivatives, s, lds, sdab, pab, force_a)
Purpose: Calculation of the two-center overlap integrals [a|b] over Cartesian Gaussian-type functions...
Definition ai_overlap.F:73
debug_os_integrals
Debugs Obara-Saika integral matrices.
Definition debug_os_integrals.F:14
debug_os_integrals::overlap_ab_test
subroutine, public overlap_ab_test(la_max, la_min, npgfa, zeta, lb_max, lb_min, npgfb, zetb, ra, rb, sab, dmax)
recursive test routines for integral (a,b)
Definition debug_os_integrals.F:135
debug_os_integrals::overlap_abc_test
subroutine, public overlap_abc_test(la_max, npgfa, zeta, la_min, lb_max, npgfb, zetb, lb_min, lc_max, npgfc, zetc, lc_min, ra, rb, rc, sabc, dmax)
recursive test routines for integral (a,b,c)
Definition debug_os_integrals.F:306
debug_os_integrals::overlap_aabb_test
subroutine, public overlap_aabb_test(la_max1, la_min1, npgfa1, zeta1, la_max2, la_min2, npgfa2, zeta2, lb_max1, lb_min1, npgfb1, zetb1, lb_max2, lb_min2, npgfb2, zetb2, ra, rb, saabb, dmax)
recursive test routines for integral (aa,bb)
Definition debug_os_integrals.F:514
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
orbital_pointers
Provides Cartesian and spherical orbital pointers and indices.
Definition orbital_pointers.F:15
orbital_pointers::ncoset
integer, dimension(:), allocatable, public ncoset
Definition orbital_pointers.F:42
orbital_pointers::coset
integer, dimension(:, :, :), allocatable, public coset
Definition orbital_pointers.F:45
orbital_pointers::indco
integer, dimension(:, :), allocatable, public indco
Definition orbital_pointers.F:43
xc
Exchange and Correlation functional calculations.
Definition xc.F:17