(git:e68414f)
Loading...
Searching...
No Matches
ai_overlap_aabb.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Calculation of the overlap integrals over Cartesian Gaussian-type
10!> functions.
11!> \par Literature
12!> S. Obara and A. Saika, J. Chem. Phys. 84, 3963 (1986)
13!> \par Parameters
14!> - ax,ay,az : Angular momentum index numbers of orbital a.
15!> - bx,by,bz : Angular momentum index numbers of orbital b.
16!> - coset : Cartesian orbital set pointer.
17!> - dab : Distance between the atomic centers a and b.
18!> - l{a,b} : Angular momentum quantum number of shell a or b.
19!> - l{a,b}_max: Maximum angular momentum quantum number of shell a or b.
20!> - l{a,b}_min: Minimum angular momentum quantum number of shell a or b.
21!> - rab : Distance vector between the atomic centers a and b.
22!> - rpgf{a,b} : Radius of the primitive Gaussian-type function a or b.
23!> - sab : Shell set of overlap integrals.
24!> - zet{a,b} : Exponents of the Gaussian-type functions a or b.
25!> - zetp : Reciprocal of the sum of the exponents of orbital a and b.
26! **************************************************************************************************
27MODULE ai_overlap_aabb
28
29 USE kinds, ONLY: dp
30 USE mathconstants, ONLY: pi
31 USE orbital_pointers, ONLY: coset,&
32 indco,&
33 ncoset
34#include "../base/base_uses.f90"
35
36 IMPLICIT NONE
37
38 PRIVATE
39
40 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'ai_overlap_aabb'
41
42! *** Public subroutines ***
43 PUBLIC :: overlap_aabb
44
45CONTAINS
46
47! **************************************************************************************************
48!> \brief Purpose: Calculation of the two-center overlap integrals [aa|bb]
49!> over Cartesian Gaussian-type functions.
50!> \param la_max_set1 ...
51!> \param la_min_set1 ...
52!> \param npgfa1 ...
53!> \param rpgfa1 ...
54!> \param zeta1 ...
55!> \param la_max_set2 ...
56!> \param la_min_set2 ...
57!> \param npgfa2 ...
58!> \param rpgfa2 ...
59!> \param zeta2 ...
60!> \param lb_max_set1 ...
61!> \param lb_min_set1 ...
62!> \param npgfb1 ...
63!> \param rpgfb1 ...
64!> \param zetb1 ...
65!> \param lb_max_set2 ...
66!> \param lb_min_set2 ...
67!> \param npgfb2 ...
68!> \param rpgfb2 ...
69!> \param zetb2 ...
70!> \param asets_equal ...
71!> \param bsets_equal ...
72!> \param rab ...
73!> \param dab ...
74!> \param saabb ...
75!> \param s ...
76!> \param lds ...
77!> \date 06.2014
78!> \author Dorothea Golze
79! **************************************************************************************************
80 SUBROUTINE overlap_aabb(la_max_set1, la_min_set1, npgfa1, rpgfa1, zeta1, &
81 la_max_set2, la_min_set2, npgfa2, rpgfa2, zeta2, &
82 lb_max_set1, lb_min_set1, npgfb1, rpgfb1, zetb1, &
83 lb_max_set2, lb_min_set2, npgfb2, rpgfb2, zetb2, &
84 asets_equal, bsets_equal, rab, dab, saabb, s, lds)
85
86 INTEGER, INTENT(IN) :: la_max_set1, la_min_set1, npgfa1
87 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: rpgfa1, zeta1
88 INTEGER, INTENT(IN) :: la_max_set2, la_min_set2, npgfa2
89 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: rpgfa2, zeta2
90 INTEGER, INTENT(IN) :: lb_max_set1, lb_min_set1, npgfb1
91 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: rpgfb1, zetb1
92 INTEGER, INTENT(IN) :: lb_max_set2, lb_min_set2, npgfb2
93 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: rpgfb2, zetb2
94 LOGICAL, INTENT(IN) :: asets_equal, bsets_equal
95 REAL(kind=dp), DIMENSION(3), INTENT(IN) :: rab
96 REAL(kind=dp), INTENT(IN) :: dab
97 REAL(kind=dp), DIMENSION(:, :, :, :), &
98 INTENT(INOUT) :: saabb
99 INTEGER, INTENT(IN) :: lds
100 REAL(kind=dp), DIMENSION(lds, lds), INTENT(INOUT) :: s
101
102 CHARACTER(len=*), PARAMETER :: routinen = 'overlap_aabb'
103
104 INTEGER :: ax, ay, az, bx, by, bz, coa, coamx, coamy, coamz, coapx, coapy, coapz, cob, &
105 cobm2x, cobm2y, cobm2z, cobmx, cobmy, cobmz, handle, i, ia, ib, ipgf, j, ja, jb, jpgf, &
106 jpgf_start, kpgf, la, la_max, la_min, la_start, lb, lb_max, lb_min, lpgf, lpgf_start, &
107 ncoa1, ncoa2, ncob1, ncob2
108 INTEGER, DIMENSION(3) :: na, naa, nb, nbb, nia, nib, nja, njb
109 REAL(kind=dp) :: f0, f1, f2, f3, f4, fax, fay, faz, zeta, &
110 zetb, zetp
111 REAL(kind=dp), DIMENSION(3) :: rap, rbp
112
113 CALL timeset(routinen, handle)
114
115! *** Loop over all pairs of primitive Gaussian-type functions ***
116
117 ncoa1 = 0
118 ncoa2 = 0
119 ncob1 = 0
120 ncob2 = 0
121
122 DO ipgf = 1, npgfa1
123
124 ncoa2 = 0
125
126 IF (asets_equal) THEN
127 jpgf_start = ipgf
128 DO i = 1, jpgf_start - 1
129 ncoa2 = ncoa2 + ncoset(la_max_set2)
130 END DO
131 ELSE
132 jpgf_start = 1
133 END IF
134
135 DO jpgf = jpgf_start, npgfa2
136
137 ncob1 = 0
138 zeta = zeta1(ipgf) + zeta2(jpgf)
139 la_max = la_max_set1 + la_max_set2
140 la_min = la_min_set1 + la_min_set2
141
142 DO kpgf = 1, npgfb1
143
144 ncob2 = 0
145
146 IF (bsets_equal) THEN
147 lpgf_start = kpgf
148 DO i = 1, lpgf_start - 1
149 ncob2 = ncob2 + ncoset(lb_max_set2)
150 END DO
151 ELSE
152 lpgf_start = 1
153 END IF
154
155 DO lpgf = lpgf_start, npgfb2
156
157 ! *** Screening ***
158 IF ((rpgfa1(ipgf) + rpgfb1(kpgf) < dab) .OR. &
159 (rpgfa2(jpgf) + rpgfb1(kpgf) < dab) .OR. &
160 (rpgfa1(ipgf) + rpgfb2(lpgf) < dab) .OR. &
161 (rpgfa2(jpgf) + rpgfb2(lpgf) < dab)) THEN
162 DO jb = ncoset(lb_min_set2 - 1) + 1, ncoset(lb_max_set2)
163 DO ib = ncoset(lb_min_set1 - 1) + 1, ncoset(lb_max_set1)
164 DO ja = ncoset(la_min_set2 - 1) + 1, ncoset(la_max_set2)
165 DO ia = ncoset(la_min_set1 - 1) + 1, ncoset(la_max_set1)
166 saabb(ncoa1 + ia, ncoa2 + ja, ncob1 + ib, ncob2 + jb) = 0._dp
167 IF (asets_equal) saabb(ncoa2 + ja, ncoa1 + ia, ncob1 + ib, ncob2 + jb) = 0._dp
168 IF (bsets_equal) saabb(ncoa1 + ia, ncoa2 + ja, ncob2 + jb, ncob1 + ib) = 0._dp
169 IF (asets_equal .AND. bsets_equal) THEN
170 saabb(ncoa2 + ja, ncoa1 + ia, ncob1 + ib, ncob2 + jb) = 0._dp
171 saabb(ncoa1 + ia, ncoa2 + ja, ncob2 + jb, ncob1 + ib) = 0._dp
172 saabb(ncoa2 + ja, ncoa1 + ia, ncob2 + jb, ncob1 + ib) = 0._dp
173 END IF
174 END DO
175 END DO
176 END DO
177 END DO
178 ncob2 = ncob2 + ncoset(lb_max_set2)
179 cycle
180 END IF
181
182 zetb = zetb1(kpgf) + zetb2(lpgf)
183 lb_max = lb_max_set1 + lb_max_set2
184 lb_min = lb_min_set1 + lb_min_set2
185
186! *** Calculate some prefactors ***
187
188 zetp = 1.0_dp/(zeta + zetb)
189
190 f0 = sqrt((pi*zetp)**3)
191 f1 = zetb*zetp
192 f2 = 0.5_dp*zetp
193
194! *** Calculate the basic two-center overlap integral [s|s] ***
195
196 s(1, 1) = f0*exp(-zeta*f1*dab*dab)
197
198! *** Recurrence steps: [s|s] -> [a|b] ***
199
200 IF (la_max > 0) THEN
201
202! *** Vertical recurrence steps: [s|s] -> [a|s] ***
203
204 rap(:) = f1*rab(:)
205
206! *** [p|s] = (Pi - Ai)*[s|s] (i = x,y,z) ***
207
208 s(2, 1) = rap(1)*s(1, 1) ! [px|s]
209 s(3, 1) = rap(2)*s(1, 1) ! [py|s]
210 s(4, 1) = rap(3)*s(1, 1) ! [pz|s]
211
212 IF (la_max > 1) THEN
213
214! *** [d|s] ***
215
216 f3 = f2*s(1, 1)
217
218 s(5, 1) = rap(1)*s(2, 1) + f3 ! [dx2|s]
219 s(6, 1) = rap(1)*s(3, 1) ! [dxy|s]
220 s(7, 1) = rap(1)*s(4, 1) ! [dxz|s]
221 s(8, 1) = rap(2)*s(3, 1) + f3 ! [dy2|s]
222 s(9, 1) = rap(2)*s(4, 1) ! [dyz|s]
223 s(10, 1) = rap(3)*s(4, 1) + f3 ! [dz2|s]
224
225 IF (la_max > 2) THEN
226
227! *** [f|s] ***
228
229 f3 = 2.0_dp*f2
230
231 s(11, 1) = rap(1)*s(5, 1) + f3*s(2, 1) ! [fx3 |s]
232 s(12, 1) = rap(1)*s(6, 1) + f2*s(3, 1) ! [fx2y|s]
233 s(13, 1) = rap(1)*s(7, 1) + f2*s(4, 1) ! [fx2z|s]
234 s(14, 1) = rap(1)*s(8, 1) ! [fxy2|s]
235 s(15, 1) = rap(1)*s(9, 1) ! [fxyz|s]
236 s(16, 1) = rap(1)*s(10, 1) ! [fxz2|s]
237 s(17, 1) = rap(2)*s(8, 1) + f3*s(3, 1) ! [fy3 |s]
238 s(18, 1) = rap(2)*s(9, 1) + f2*s(4, 1) ! [fy2z|s]
239 s(19, 1) = rap(2)*s(10, 1) ! [fyz2|s]
240 s(20, 1) = rap(3)*s(10, 1) + f3*s(4, 1) ! [fz3 |s]
241
242 IF (la_max > 3) THEN
243
244! *** [g|s] ***
245
246 f4 = 3.0_dp*f2
247
248 s(21, 1) = rap(1)*s(11, 1) + f4*s(5, 1) ! [gx4 |s]
249 s(22, 1) = rap(1)*s(12, 1) + f3*s(6, 1) ! [gx3y |s]
250 s(23, 1) = rap(1)*s(13, 1) + f3*s(7, 1) ! [gx3z |s]
251 s(24, 1) = rap(1)*s(14, 1) + f2*s(8, 1) ! [gx2y2|s]
252 s(25, 1) = rap(1)*s(15, 1) + f2*s(9, 1) ! [gx2yz|s]
253 s(26, 1) = rap(1)*s(16, 1) + f2*s(10, 1) ! [gx2z2|s]
254 s(27, 1) = rap(1)*s(17, 1) ! [gxy3 |s]
255 s(28, 1) = rap(1)*s(18, 1) ! [gxy2z|s]
256 s(29, 1) = rap(1)*s(19, 1) ! [gxyz2|s]
257 s(30, 1) = rap(1)*s(20, 1) ! [gxz3 |s]
258 s(31, 1) = rap(2)*s(17, 1) + f4*s(8, 1) ! [gy4 |s]
259 s(32, 1) = rap(2)*s(18, 1) + f3*s(9, 1) ! [gy3z |s]
260 s(33, 1) = rap(2)*s(19, 1) + f2*s(10, 1) ! [gy2z2|s]
261 s(34, 1) = rap(2)*s(20, 1) ! [gyz3 |s]
262 s(35, 1) = rap(3)*s(20, 1) + f4*s(10, 1) ! [gz4 |s]
263
264! *** [a|s] = (Pi - Ai)*[a-1i|s] + f2*Ni(a-1i)*[a-2i|s] ***
265
266 DO la = 5, la_max
267
268! *** Increase the angular momentum component z of a ***
269
270 s(coset(0, 0, la), 1) = &
271 rap(3)*s(coset(0, 0, la - 1), 1) + &
272 f2*real(la - 1, dp)*s(coset(0, 0, la - 2), 1)
273
274! *** Increase the angular momentum component y of a ***
275
276 az = la - 1
277 s(coset(0, 1, az), 1) = rap(2)*s(coset(0, 0, az), 1)
278 DO ay = 2, la
279 az = la - ay
280 s(coset(0, ay, az), 1) = &
281 rap(2)*s(coset(0, ay - 1, az), 1) + &
282 f2*real(ay - 1, dp)*s(coset(0, ay - 2, az), 1)
283 END DO
284
285! *** Increase the angular momentum component x of a ***
286
287 DO ay = 0, la - 1
288 az = la - 1 - ay
289 s(coset(1, ay, az), 1) = rap(1)*s(coset(0, ay, az), 1)
290 END DO
291 DO ax = 2, la
292 f3 = f2*real(ax - 1, dp)
293 DO ay = 0, la - ax
294 az = la - ax - ay
295 s(coset(ax, ay, az), 1) = &
296 rap(1)*s(coset(ax - 1, ay, az), 1) + &
297 f3*s(coset(ax - 2, ay, az), 1)
298 END DO
299 END DO
300
301 END DO
302
303 END IF
304
305 END IF
306
307 END IF
308
309! *** Recurrence steps: [a|s] -> [a|b] ***
310
311 IF (lb_max > 0) THEN
312
313 DO j = 2, ncoset(lb_max)
314 DO i = 1, ncoset(la_min)
315 s(i, j) = 0.0_dp
316 END DO
317 END DO
318
319! *** Horizontal recurrence steps ***
320
321 rbp(:) = rap(:) - rab(:)
322
323! *** [a|p] = [a+1i|s] - (Bi - Ai)*[a|s] ***
324
325 IF (lb_max == 1) THEN
326 la_start = la_min
327 ELSE
328 la_start = max(0, la_min - 1)
329 END IF
330
331 DO la = la_start, la_max - 1
332 DO ax = 0, la
333 DO ay = 0, la - ax
334 az = la - ax - ay
335 coa = coset(ax, ay, az)
336 coapx = coset(ax + 1, ay, az)
337 coapy = coset(ax, ay + 1, az)
338 coapz = coset(ax, ay, az + 1)
339 s(coa, 2) = s(coapx, 1) - rab(1)*s(coa, 1)
340 s(coa, 3) = s(coapy, 1) - rab(2)*s(coa, 1)
341 s(coa, 4) = s(coapz, 1) - rab(3)*s(coa, 1)
342 END DO
343 END DO
344 END DO
345
346! *** Vertical recurrence step ***
347
348! *** [a|p] = (Pi - Bi)*[a|s] + f2*Ni(a)*[a-1i|s] ***
349
350 DO ax = 0, la_max
351 fax = f2*real(ax, dp)
352 DO ay = 0, la_max - ax
353 fay = f2*real(ay, dp)
354 az = la_max - ax - ay
355 faz = f2*real(az, dp)
356 coa = coset(ax, ay, az)
357 coamx = coset(ax - 1, ay, az)
358 coamy = coset(ax, ay - 1, az)
359 coamz = coset(ax, ay, az - 1)
360 s(coa, 2) = rbp(1)*s(coa, 1) + fax*s(coamx, 1)
361 s(coa, 3) = rbp(2)*s(coa, 1) + fay*s(coamy, 1)
362 s(coa, 4) = rbp(3)*s(coa, 1) + faz*s(coamz, 1)
363 END DO
364 END DO
365
366! *** Recurrence steps: [a|p] -> [a|b] ***
367
368 DO lb = 2, lb_max
369
370! *** Horizontal recurrence steps ***
371
372! *** [a|b] = [a+1i|b-1i] - (Bi - Ai)*[a|b-1i] ***
373
374 IF (lb == lb_max) THEN
375 la_start = la_min
376 ELSE
377 la_start = max(0, la_min - 1)
378 END IF
379
380 DO la = la_start, la_max - 1
381 DO ax = 0, la
382 DO ay = 0, la - ax
383 az = la - ax - ay
384 coa = coset(ax, ay, az)
385 coapx = coset(ax + 1, ay, az)
386 coapy = coset(ax, ay + 1, az)
387 coapz = coset(ax, ay, az + 1)
388
389! *** Shift of angular momentum component z from a to b ***
390
391 cob = coset(0, 0, lb)
392 cobmz = coset(0, 0, lb - 1)
393 s(coa, cob) = s(coapz, cobmz) - rab(3)*s(coa, cobmz)
394
395! *** Shift of angular momentum component y from a to b ***
396
397 DO by = 1, lb
398 bz = lb - by
399 cob = coset(0, by, bz)
400 cobmy = coset(0, by - 1, bz)
401 s(coa, cob) = s(coapy, cobmy) - rab(2)*s(coa, cobmy)
402 END DO
403
404! *** Shift of angular momentum component x from a to b ***
405
406 DO bx = 1, lb
407 DO by = 0, lb - bx
408 bz = lb - bx - by
409 cob = coset(bx, by, bz)
410 cobmx = coset(bx - 1, by, bz)
411 s(coa, cob) = s(coapx, cobmx) - rab(1)*s(coa, cobmx)
412 END DO
413 END DO
414
415 END DO
416 END DO
417 END DO
418
419! *** Vertical recurrence step ***
420
421! *** [a|b] = (Pi - Bi)*[a|b-1i] + f2*Ni(a)*[a-1i|b-1i] + ***
422! *** f2*Ni(b-1i)*[a|b-2i] ***
423
424 DO ax = 0, la_max
425 fax = f2*real(ax, dp)
426 DO ay = 0, la_max - ax
427 fay = f2*real(ay, dp)
428 az = la_max - ax - ay
429 faz = f2*real(az, dp)
430 coa = coset(ax, ay, az)
431 coamx = coset(ax - 1, ay, az)
432 coamy = coset(ax, ay - 1, az)
433 coamz = coset(ax, ay, az - 1)
434
435! *** Increase the angular momentum component z of b ***
436
437 f3 = f2*real(lb - 1, dp)
438 cob = coset(0, 0, lb)
439 cobmz = coset(0, 0, lb - 1)
440 cobm2z = coset(0, 0, lb - 2)
441 s(coa, cob) = rbp(3)*s(coa, cobmz) + &
442 faz*s(coamz, cobmz) + &
443 f3*s(coa, cobm2z)
444
445! *** Increase the angular momentum component y of b ***
446
447 bz = lb - 1
448 cob = coset(0, 1, bz)
449 cobmy = coset(0, 0, bz)
450 s(coa, cob) = rbp(2)*s(coa, cobmy) + &
451 fay*s(coamy, cobmy)
452 DO by = 2, lb
453 bz = lb - by
454 f3 = f2*real(by - 1, dp)
455 cob = coset(0, by, bz)
456 cobmy = coset(0, by - 1, bz)
457 cobm2y = coset(0, by - 2, bz)
458 s(coa, cob) = rbp(2)*s(coa, cobmy) + &
459 fay*s(coamy, cobmy) + &
460 f3*s(coa, cobm2y)
461 END DO
462
463! *** Increase the angular momentum component x of b ***
464
465 DO by = 0, lb - 1
466 bz = lb - 1 - by
467 cob = coset(1, by, bz)
468 cobmx = coset(0, by, bz)
469 s(coa, cob) = rbp(1)*s(coa, cobmx) + &
470 fax*s(coamx, cobmx)
471 END DO
472 DO bx = 2, lb
473 f3 = f2*real(bx - 1, dp)
474 DO by = 0, lb - bx
475 bz = lb - bx - by
476 cob = coset(bx, by, bz)
477 cobmx = coset(bx - 1, by, bz)
478 cobm2x = coset(bx - 2, by, bz)
479 s(coa, cob) = rbp(1)*s(coa, cobmx) + &
480 fax*s(coamx, cobmx) + &
481 f3*s(coa, cobm2x)
482 END DO
483 END DO
484
485 END DO
486 END DO
487
488 END DO
489
490 END IF
491
492 ELSE
493
494 IF (lb_max > 0) THEN
495
496! *** Vertical recurrence steps: [s|s] -> [s|b] ***
497
498 rbp(:) = (f1 - 1.0_dp)*rab(:)
499
500! *** [s|p] = (Pi - Bi)*[s|s] ***
501
502 s(1, 2) = rbp(1)*s(1, 1) ! [s|px]
503 s(1, 3) = rbp(2)*s(1, 1) ! [s|py]
504 s(1, 4) = rbp(3)*s(1, 1) ! [s|pz]
505
506 IF (lb_max > 1) THEN
507
508! *** [s|d] ***
509
510 f3 = f2*s(1, 1)
511
512 s(1, 5) = rbp(1)*s(1, 2) + f3 ! [s|dx2]
513 s(1, 6) = rbp(1)*s(1, 3) ! [s|dxy]
514 s(1, 7) = rbp(1)*s(1, 4) ! [s|dxz]
515 s(1, 8) = rbp(2)*s(1, 3) + f3 ! [s|dy2]
516 s(1, 9) = rbp(2)*s(1, 4) ! [s|dyz]
517 s(1, 10) = rbp(3)*s(1, 4) + f3 ! [s|dz2]
518
519! *** [s|b] = (Pi - Bi)*[s|b-1i] + f2*Ni(b-1i)*[s|b-2i] ***
520
521 DO lb = 3, lb_max
522
523! *** Increase the angular momentum component z of b ***
524
525 s(1, coset(0, 0, lb)) = &
526 rbp(3)*s(1, coset(0, 0, lb - 1)) + &
527 f2*real(lb - 1, dp)*s(1, coset(0, 0, lb - 2))
528
529! *** Increase the angular momentum component y of b ***
530
531 bz = lb - 1
532 s(1, coset(0, 1, bz)) = rbp(2)*s(1, coset(0, 0, bz))
533 DO by = 2, lb
534 bz = lb - by
535 s(1, coset(0, by, bz)) = &
536 rbp(2)*s(1, coset(0, by - 1, bz)) + &
537 f2*real(by - 1, dp)*s(1, coset(0, by - 2, bz))
538 END DO
539
540! *** Increase the angular momentum component x of b ***
541
542 DO by = 0, lb - 1
543 bz = lb - 1 - by
544 s(1, coset(1, by, bz)) = rbp(1)*s(1, coset(0, by, bz))
545 END DO
546 DO bx = 2, lb
547 f3 = f2*real(bx - 1, dp)
548 DO by = 0, lb - bx
549 bz = lb - bx - by
550 s(1, coset(bx, by, bz)) = &
551 rbp(1)*s(1, coset(bx - 1, by, bz)) + &
552 f3*s(1, coset(bx - 2, by, bz))
553 END DO
554 END DO
555
556 END DO
557
558 END IF
559
560 END IF
561
562 END IF
563
564! *** Store the primitive overlap integrals ***
565 DO jb = ncoset(lb_min_set2 - 1) + 1, ncoset(lb_max_set2)
566 njb(1:3) = indco(1:3, jb)
567 DO ib = ncoset(lb_min_set1 - 1) + 1, ncoset(lb_max_set1)
568 nib(1:3) = indco(1:3, ib)
569 nbb(1:3) = nib + njb
570 DO ja = ncoset(la_min_set2 - 1) + 1, ncoset(la_max_set2)
571 nja(1:3) = indco(1:3, ja)
572 DO ia = ncoset(la_min_set1 - 1) + 1, ncoset(la_max_set1)
573 nia(1:3) = indco(1:3, ia)
574 naa(1:3) = nia + nja
575 ! now loop over all elements of s
576 DO j = ncoset(lb_min - 1) + 1, ncoset(lb_max)
577 nb(1:3) = indco(1:3, j)
578 DO i = ncoset(la_min - 1) + 1, ncoset(la_max)
579 na(1:3) = indco(1:3, i)
580 IF (all(na == naa) .AND. all(nb == nbb)) THEN
581 saabb(ncoa1 + ia, ncoa2 + ja, ncob1 + ib, ncob2 + jb) = s(i, j)
582 IF (asets_equal) saabb(ncoa2 + ja, ncoa1 + ia, ncob1 + ib, ncob2 + jb) = s(i, j)
583 IF (bsets_equal) saabb(ncoa1 + ia, ncoa2 + ja, ncob2 + jb, ncob1 + ib) = s(i, j)
584 IF (asets_equal .AND. bsets_equal) THEN
585 saabb(ncoa2 + ja, ncoa1 + ia, ncob1 + ib, ncob2 + jb) = s(i, j)
586 saabb(ncoa1 + ia, ncoa2 + ja, ncob2 + jb, ncob1 + ib) = s(i, j)
587 saabb(ncoa2 + ja, ncoa1 + ia, ncob2 + jb, ncob1 + ib) = s(i, j)
588 END IF
589 END IF
590 END DO
591 END DO
592 END DO
593 END DO
594 END DO
595 END DO
596
597 ncob2 = ncob2 + ncoset(lb_max_set2)
598
599 END DO
600
601 ncob1 = ncob1 + ncoset(lb_max_set1)
602
603 END DO
604
605 ncoa2 = ncoa2 + ncoset(la_max_set2)
606
607 END DO
608
609 ncoa1 = ncoa1 + ncoset(la_max_set1)
610
611 END DO
612
613 CALL timestop(handle)
614
615 END SUBROUTINE overlap_aabb
616
617END MODULE ai_overlap_aabb
ai_overlap_aabb
Calculation of the overlap integrals over Cartesian Gaussian-type functions.
Definition ai_overlap_aabb.F:27
ai_overlap_aabb::overlap_aabb
subroutine, public overlap_aabb(la_max_set1, la_min_set1, npgfa1, rpgfa1, zeta1, la_max_set2, la_min_set2, npgfa2, rpgfa2, zeta2, lb_max_set1, lb_min_set1, npgfb1, rpgfb1, zetb1, lb_max_set2, lb_min_set2, npgfb2, rpgfb2, zetb2, asets_equal, bsets_equal, rab, dab, saabb, s, lds)
Purpose: Calculation of the two-center overlap integrals [aa|bb] over Cartesian Gaussian-type functio...
Definition ai_overlap_aabb.F:85
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::pi
real(kind=dp), parameter, public pi
Definition mathconstants.F:110
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