(git:0761fd1)
Loading...
Searching...
No Matches
ai_os_rr.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8MODULE ai_os_rr
9
10 USE gamma, ONLY: fgamma => fgamma_0
11 USE kinds, ONLY: dp
12 USE orbital_pointers, ONLY: coset, indco, ncoset
13#include "../base/base_uses.f90"
14
15 IMPLICIT NONE
16 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'ai_os_rr'
17 PRIVATE
18
19 ! *** Public subroutines ***
20
21 PUBLIC :: os_rr_ovlp, os_rr_ovlp3, os_rr_coul, os_rr_radial, os_rr_coulomb3
22
23 INTERFACE os_rr_ovlp
24 MODULE PROCEDURE os_rr_ovlp_real, os_rr_ovlp_complex
25 END INTERFACE os_rr_ovlp
26
27 INTERFACE os_rr_ovlp3
28 MODULE PROCEDURE os_rr_ovlp3_real, os_rr_ovlp3_complex
29 END INTERFACE os_rr_ovlp3
30
31 INTERFACE os_rr_radial
32 MODULE PROCEDURE os_rr_radial_abm
33 END INTERFACE os_rr_radial
34
35CONTAINS
36
37! **************************************************************************************************
38!> \brief Calculation of the basic Obara-Saika recurrence relation
39!> for real or complex product-center displacements. Each Cartesian
40!> table is normalized to a unit zeroth-degree integral.
41!> \param rap product-center displacement relative to center A
42!> \param la_max ...
43!> \param rbp product-center displacement relative to center B
44!> \param lb_max ...
45!> \param zet real sum of the Gaussian exponents
46!> \param ldrr ...
47!> \param rr ...
48!> \date 02.03.2009
49!> \author VW
50!> \version 1.0
51! **************************************************************************************************
52 SUBROUTINE os_rr_ovlp_real (rap, la_max, rbp, lb_max, zet, ldrr, rr)
53 REAL (dp), DIMENSION(3), INTENT(IN) :: rap
54 INTEGER, INTENT(IN) :: la_max
55 REAL (dp), DIMENSION(3), INTENT(IN) :: rbp
56 INTEGER, INTENT(IN) :: lb_max
57 REAL(dp), INTENT(IN) :: zet
58 INTEGER, INTENT(IN) :: ldrr
59 REAL (dp), DIMENSION(0:ldrr - 1, 0:ldrr - 1, 3) :: rr
60
61 INTEGER :: idir
62
63 DO idir = 1, 3
64 CALL os_rr_ovlp_1d_real (rap(idir), la_max, rbp(idir), lb_max, zet, rr(:, :, idir))
65 END DO
66 END SUBROUTINE os_rr_ovlp_real
67
68! **************************************************************************************************
69!> \brief Normalized one-dimensional two-center Gaussian moments.
70!> \param rap product-center displacement relative to A
71!> \param la_max highest A degree
72!> \param rbp product-center displacement relative to B
73!> \param lb_max highest B degree
74!> \param zet sum of Gaussian exponents
75!> \param rr normalized moments with zeroth-degree value one
76! **************************************************************************************************
77 PURE SUBROUTINE os_rr_ovlp_1d_real (rap, la_max, rbp, lb_max, zet, rr)
78 REAL (dp), INTENT(IN) :: rap, rbp
79 INTEGER, INTENT(IN) :: la_max, lb_max
80 REAL(dp), INTENT(IN) :: zet
81 REAL (dp), INTENT(INOUT) :: rr(0:, 0:)
82
83 INTEGER :: la, lam1, lap1, lb, lbm1, lbp1
84 REAL(dp) :: g
85
86 g = 0.5_dp/zet
87 rr(0, 0) = 1.0_dp
88 ! recursion along la for lb=0
89 IF (la_max > 0) THEN
90 rr(1, 0) = rap
91 DO la = 1, la_max - 1
92 lap1 = la + 1
93 lam1 = la - 1
94 rr(lap1, 0) = real(la, dp)*g*rr(lam1, 0) + rap*rr(la, 0)
95 END DO
96 END IF
97 ! recursion along lb for all la
98 IF (lb_max > 0) THEN
99 rr(0, 1) = rbp
100 DO la = 1, la_max
101 lam1 = la - 1
102 rr(la, 1) = real(la, dp)*g*rr(lam1, 0) + rbp*rr(la, 0)
103 END DO
104 DO lb = 1, lb_max - 1
105 lbp1 = lb + 1
106 lbm1 = lb - 1
107 rr(0, lbp1) = real(lb, dp)*g*rr(0, lbm1) + rbp*rr(0, lb)
108 DO la = 1, la_max
109 lam1 = la - 1
110 rr(la, lbp1) = g*(real(la, dp)*rr(lam1, lb) + real(lb, dp)*rr(la, lbm1)) + &
111 rbp*rr(la, lb)
112 END DO
113 END DO
114 END IF
115 END SUBROUTINE os_rr_ovlp_1d_real
116
117! **************************************************************************************************
118!> \brief Normalized separable three-center Gaussian overlap tables.
119!> The C recurrence uses the directly supplied G-C displacement.
120!> \param rap product-center displacement G-A
121!> \param la_max highest A degree
122!> \param rbp product-center displacement G-B
123!> \param lb_max highest B degree
124!> \param rcp product-center displacement G-C
125!> \param lc_max highest C degree
126!> \param zet sum of all three Gaussian exponents
127!> \param rr normalized tables indexed by A/B/C degrees and Cartesian axis
128! **************************************************************************************************
129 SUBROUTINE os_rr_ovlp3_real (rap, la_max, rbp, lb_max, rcp, lc_max, zet, rr)
130 REAL (dp), INTENT(IN) :: rap(3), rbp(3), rcp(3)
131 INTEGER, INTENT(IN) :: la_max, lb_max, lc_max
132 REAL(dp), INTENT(IN) :: zet
133 REAL (dp), INTENT(INOUT) :: rr(0:, 0:, 0:, :)
134
135 INTEGER :: a, b, c, idir
136 REAL(dp) :: g
137
138 g = 0.5_dp/zet
139 DO idir = 1, 3
140 CALL os_rr_ovlp_1d_real (rap(idir), la_max, rbp(idir), lb_max, zet, rr(:, :, 0, idir))
141 DO c = 0, lc_max - 1
142 DO b = 0, lb_max
143 DO a = 0, la_max
144 rr(a, b, c + 1, idir) = rcp(idir)*rr(a, b, c, idir) + &
145 g*(real(a, dp)*rr(max(0, a - 1), b, c, idir) + &
146 REAL(b, dp)*rr(a, MAX(0, b - 1), c, idir) + &
147 REAL(c, dp)*rr(a, b, MAX(0, c - 1), idir))
148 END DO
149 END DO
150 END DO
151 END DO
152 END SUBROUTINE os_rr_ovlp3_real
153! **************************************************************************************************
154!> \brief Calculation of the basic Obara-Saika recurrence relation
155!> for real or complex product-center displacements. Each Cartesian
156!> table is normalized to a unit zeroth-degree integral.
157!> \param rap product-center displacement relative to center A
158!> \param la_max ...
159!> \param rbp product-center displacement relative to center B
160!> \param lb_max ...
161!> \param zet real sum of the Gaussian exponents
162!> \param ldrr ...
163!> \param rr ...
164!> \date 02.03.2009
165!> \author VW
166!> \version 1.0
167! **************************************************************************************************
168 SUBROUTINE os_rr_ovlp_complex (rap, la_max, rbp, lb_max, zet, ldrr, rr)
169 COMPLEX (dp), DIMENSION(3), INTENT(IN) :: rap
170 INTEGER, INTENT(IN) :: la_max
171 COMPLEX (dp), DIMENSION(3), INTENT(IN) :: rbp
172 INTEGER, INTENT(IN) :: lb_max
173 REAL(dp), INTENT(IN) :: zet
174 INTEGER, INTENT(IN) :: ldrr
175 COMPLEX (dp), DIMENSION(0:ldrr - 1, 0:ldrr - 1, 3) :: rr
176
177 INTEGER :: idir
178
179 DO idir = 1, 3
180 CALL os_rr_ovlp_1d_complex (rap(idir), la_max, rbp(idir), lb_max, zet, rr(:, :, idir))
181 END DO
182 END SUBROUTINE os_rr_ovlp_complex
183
184! **************************************************************************************************
185!> \brief Normalized one-dimensional two-center Gaussian moments.
186!> \param rap product-center displacement relative to A
187!> \param la_max highest A degree
188!> \param rbp product-center displacement relative to B
189!> \param lb_max highest B degree
190!> \param zet sum of Gaussian exponents
191!> \param rr normalized moments with zeroth-degree value one
192! **************************************************************************************************
193 PURE SUBROUTINE os_rr_ovlp_1d_complex (rap, la_max, rbp, lb_max, zet, rr)
194 COMPLEX (dp), INTENT(IN) :: rap, rbp
195 INTEGER, INTENT(IN) :: la_max, lb_max
196 REAL(dp), INTENT(IN) :: zet
197 COMPLEX (dp), INTENT(INOUT) :: rr(0:, 0:)
198
199 INTEGER :: la, lam1, lap1, lb, lbm1, lbp1
200 REAL(dp) :: g
201
202 g = 0.5_dp/zet
203 rr(0, 0) = 1.0_dp
204 ! recursion along la for lb=0
205 IF (la_max > 0) THEN
206 rr(1, 0) = rap
207 DO la = 1, la_max - 1
208 lap1 = la + 1
209 lam1 = la - 1
210 rr(lap1, 0) = real(la, dp)*g*rr(lam1, 0) + rap*rr(la, 0)
211 END DO
212 END IF
213 ! recursion along lb for all la
214 IF (lb_max > 0) THEN
215 rr(0, 1) = rbp
216 DO la = 1, la_max
217 lam1 = la - 1
218 rr(la, 1) = real(la, dp)*g*rr(lam1, 0) + rbp*rr(la, 0)
219 END DO
220 DO lb = 1, lb_max - 1
221 lbp1 = lb + 1
222 lbm1 = lb - 1
223 rr(0, lbp1) = real(lb, dp)*g*rr(0, lbm1) + rbp*rr(0, lb)
224 DO la = 1, la_max
225 lam1 = la - 1
226 rr(la, lbp1) = g*(real(la, dp)*rr(lam1, lb) + real(lb, dp)*rr(la, lbm1)) + &
227 rbp*rr(la, lb)
228 END DO
229 END DO
230 END IF
231 END SUBROUTINE os_rr_ovlp_1d_complex
232
233! **************************************************************************************************
234!> \brief Normalized separable three-center Gaussian overlap tables.
235!> The C recurrence uses the directly supplied G-C displacement.
236!> \param rap product-center displacement G-A
237!> \param la_max highest A degree
238!> \param rbp product-center displacement G-B
239!> \param lb_max highest B degree
240!> \param rcp product-center displacement G-C
241!> \param lc_max highest C degree
242!> \param zet sum of all three Gaussian exponents
243!> \param rr normalized tables indexed by A/B/C degrees and Cartesian axis
244! **************************************************************************************************
245 SUBROUTINE os_rr_ovlp3_complex (rap, la_max, rbp, lb_max, rcp, lc_max, zet, rr)
246 COMPLEX (dp), INTENT(IN) :: rap(3), rbp(3), rcp(3)
247 INTEGER, INTENT(IN) :: la_max, lb_max, lc_max
248 REAL(dp), INTENT(IN) :: zet
249 COMPLEX (dp), INTENT(INOUT) :: rr(0:, 0:, 0:, :)
250
251 INTEGER :: a, b, c, idir
252 REAL(dp) :: g
253
254 g = 0.5_dp/zet
255 DO idir = 1, 3
256 CALL os_rr_ovlp_1d_complex (rap(idir), la_max, rbp(idir), lb_max, zet, rr(:, :, 0, idir))
257 DO c = 0, lc_max - 1
258 DO b = 0, lb_max
259 DO a = 0, la_max
260 rr(a, b, c + 1, idir) = rcp(idir)*rr(a, b, c, idir) + &
261 g*(real(a, dp)*rr(max(0, a - 1), b, c, idir) + &
262 REAL(b, dp)*rr(a, MAX(0, b - 1), c, idir) + &
263 REAL(c, dp)*rr(a, b, MAX(0, c - 1), idir))
264 END DO
265 END DO
266 END DO
267 END DO
268 END SUBROUTINE os_rr_ovlp3_complex
269
270! **************************************************************************************************
271!> \brief Propagate caller-supplied radial seeds through the A/B angular recurrence.
272!> Both supported layouts use the same algebra and require no transpose.
273!> \param rap displacement P-A
274!> \param la_max highest A shell
275!> \param rbp displacement P-B
276!> \param lb_max highest B shell
277!> \param rpw displacement multiplying the next radial order
278!> \param zetp reciprocal exponent sum on A/B
279!> \param eta reduction factor of the next radial order (one for point Coulomb)
280!> \param mmax highest supplied seed order, including subsequent C recurrence
281!> \param rr caller-initialized s/s seeds and resulting angular table
282! **************************************************************************************************
283 SUBROUTINE os_rr_radial_abm (rap, la_max, rbp, lb_max, rpw, zetp, eta, mmax, rr)
284 REAL(dp), INTENT(IN) :: rap(3), rbp(3), rpw(3), zetp, eta
285 INTEGER, INTENT(IN) :: la_max, lb_max, mmax
286 REAL(dp), INTENT(INOUT) :: rr(:, :, 0:)
287
288 INTEGER :: axis, degree, ia, ia1, ia2, ib, ib1, ib2, ico, iinner, iouter, la, lb, m, power(3)
289 INTEGER :: lowered(3, ncoset(MAX(la_max, lb_max)))
290 REAL(dp) :: g
291
292 g = 0.5_dp*zetp
293 DO ico = 1, ncoset(max(la_max, lb_max))
294 power = indco(:, ico)
295 lowered(:, ico) = [coset(max(0, power(1) - 1), power(2), power(3)), &
296 coset(power(1), max(0, power(2) - 1), power(3)), &
297 coset(power(1), power(2), max(0, power(3) - 1))]
298 END DO
299 DO la = 1, la_max
300 DO iouter = 0, mmax - la
301 DO iinner = ncoset(la - 1) + 1, ncoset(la)
302 m = iouter
303 ia = iinner
304 axis = 3
305 IF (indco(3, ia) == 0) axis = 2
306 IF (indco(axis, ia) == 0) axis = 1
307 ia1 = lowered(axis, ia)
308 ia2 = lowered(axis, ia1)
309 degree = indco(axis, ia) - 1
310 rr(ia, 1, m) = rap(axis)*rr(ia1, 1, m) + &
311 rpw(axis)*rr(ia1, 1, m + 1)
312 IF (degree > 0) THEN
313 rr(ia, 1, m) = rr(ia, 1, m) + g*real(degree, dp)* &
314 (rr(ia2, 1, m) - eta*rr(ia2, 1, m + 1))
315 END IF
316 END DO
317 END DO
318 END DO
319 DO lb = 1, lb_max
320 DO ib = ncoset(lb - 1) + 1, ncoset(lb)
321 axis = 3
322 IF (indco(3, ib) == 0) axis = 2
323 IF (indco(axis, ib) == 0) axis = 1
324 ib1 = lowered(axis, ib)
325 ib2 = lowered(axis, ib1)
326 degree = indco(axis, ib) - 1
327 DO la = 0, la_max
328 DO iouter = 0, mmax - la - lb
329 DO iinner = ncoset(la - 1) + 1, ncoset(la)
330 m = iouter
331 ia = iinner
332 ia1 = lowered(axis, ia)
333 rr(ia, ib, m) = rbp(axis)*rr(ia, ib1, m) + &
334 rpw(axis)*rr(ia, ib1, m + 1)
335 IF (degree > 0) THEN
336 rr(ia, ib, m) = rr(ia, ib, m) + g*real(degree, dp)* &
337 (rr(ia, ib2, m) - eta*rr(ia, ib2, m + 1))
338 END IF
339 IF (indco(axis, ia) > 0) THEN
340 rr(ia, ib, m) = rr(ia, ib, m) + g*real(indco(axis, ia), dp)* &
341 (rr(ia1, ib1, m) - eta*rr(ia1, ib1, m + 1))
342 END IF
343 END DO
344 END DO
345 END DO
346 END DO
347 END DO
348 END SUBROUTINE os_rr_radial_abm
349! **************************************************************************************************
350!> \brief Propagate caller-supplied radial seeds through the A/B angular recurrence.
351!> Both supported layouts use the same algebra and require no transpose.
352!> \param rap displacement P-A
353!> \param la_max highest A shell
354!> \param rbp displacement P-B
355!> \param lb_max highest B shell
356!> \param rpw displacement multiplying the next radial order
357!> \param zetp reciprocal exponent sum on A/B
358!> \param eta reduction factor of the next radial order (one for point Coulomb)
359!> \param mmax highest supplied seed order, including subsequent C recurrence
360!> \param rr caller-initialized s/s seeds and resulting angular table
361! **************************************************************************************************
362 SUBROUTINE os_rr_radial_mab (rap, la_max, rbp, lb_max, rpw, zetp, eta, mmax, rr)
363 REAL(dp), INTENT(IN) :: rap(3), rbp(3), rpw(3), zetp, eta
364 INTEGER, INTENT(IN) :: la_max, lb_max, mmax
365 REAL(dp), INTENT(INOUT) :: rr(0:, :, :)
366
367 INTEGER :: axis, degree, ia, ia1, ia2, ib, ib1, ib2, ico, iinner, iouter, la, lb, m, power(3)
368 INTEGER :: lowered(3, ncoset(MAX(la_max, lb_max)))
369 REAL(dp) :: g
370
371 g = 0.5_dp*zetp
372 DO ico = 1, ncoset(max(la_max, lb_max))
373 power = indco(:, ico)
374 lowered(:, ico) = [coset(max(0, power(1) - 1), power(2), power(3)), &
375 coset(power(1), max(0, power(2) - 1), power(3)), &
376 coset(power(1), power(2), max(0, power(3) - 1))]
377 END DO
378 DO la = 1, la_max
379 DO iouter = ncoset(la - 1) + 1, ncoset(la)
380 DO iinner = 0, mmax - la
381 ia = iouter
382 m = iinner
383 axis = 3
384 IF (indco(3, ia) == 0) axis = 2
385 IF (indco(axis, ia) == 0) axis = 1
386 ia1 = lowered(axis, ia)
387 ia2 = lowered(axis, ia1)
388 degree = indco(axis, ia) - 1
389 rr(m, ia, 1) = rap(axis)*rr(m, ia1, 1) + &
390 rpw(axis)*rr(m + 1, ia1, 1)
391 IF (degree > 0) THEN
392 rr(m, ia, 1) = rr(m, ia, 1) + g*real(degree, dp)* &
393 (rr(m, ia2, 1) - eta*rr(m + 1, ia2, 1))
394 END IF
395 END DO
396 END DO
397 END DO
398 DO lb = 1, lb_max
399 DO ib = ncoset(lb - 1) + 1, ncoset(lb)
400 axis = 3
401 IF (indco(3, ib) == 0) axis = 2
402 IF (indco(axis, ib) == 0) axis = 1
403 ib1 = lowered(axis, ib)
404 ib2 = lowered(axis, ib1)
405 degree = indco(axis, ib) - 1
406 DO la = 0, la_max
407 DO iouter = ncoset(la - 1) + 1, ncoset(la)
408 DO iinner = 0, mmax - la - lb
409 ia = iouter
410 m = iinner
411 ia1 = lowered(axis, ia)
412 rr(m, ia, ib) = rbp(axis)*rr(m, ia, ib1) + &
413 rpw(axis)*rr(m + 1, ia, ib1)
414 IF (degree > 0) THEN
415 rr(m, ia, ib) = rr(m, ia, ib) + g*real(degree, dp)* &
416 (rr(m, ia, ib2) - eta*rr(m + 1, ia, ib2))
417 END IF
418 IF (indco(axis, ia) > 0) THEN
419 rr(m, ia, ib) = rr(m, ia, ib) + g*real(indco(axis, ia), dp)* &
420 (rr(m, ia1, ib1) - eta*rr(m + 1, ia1, ib1))
421 END IF
422 END DO
423 END DO
424 END DO
425 END DO
426 END DO
427 END SUBROUTINE os_rr_radial_mab
428
429! **************************************************************************************************
430!> \brief Angular three-center Coulomb recurrence from caller-supplied s/s/s seeds.
431!> \param rap displacement P-A
432!> \param la_max highest A shell
433!> \param rbp displacement P-B
434!> \param lb_max highest B shell
435!> \param rcp displacement P-C
436!> \param lc_max highest C shell
437!> \param zetp reciprocal exponent sum on A/B
438!> \param zetq reciprocal exponent on C
439!> \param zetw reciprocal exponent sum on A/B/C
440!> \param rr caller-initialized seeds and angular table, with radial order last
441! **************************************************************************************************
442 SUBROUTINE os_rr_coulomb3(rap, la_max, rbp, lb_max, rcp, lc_max, zetp, zetq, zetw, rr)
443 REAL(dp), INTENT(IN) :: rap(3), rbp(3), rcp(3), zetp, zetq, zetw
444 INTEGER, INTENT(IN) :: la_max, lb_max, lc_max
445 REAL(dp), INTENT(INOUT) :: rr(:, :, :, 0:)
446
447 INTEGER :: axis, ia, ia1, ib, ib1, ic, ic1, ic2, ico, la, lb, lc, m, mmax, power(3)
448 INTEGER :: lowered(3, ncoset(max(la_max, lb_max, lc_max)))
449 REAL(dp) :: eta, eta_c, g, g_c, rcw(3), rpw(3)
450
451 mmax = la_max + lb_max + lc_max
452 eta = zetw/zetq
453 eta_c = zetw/zetp
454 rpw = -eta*rcp
455 ! Direct evaluation avoids subtracting P-C and its nearly equal Coulomb shift.
456 rcw = eta_c*rcp
457 CALL os_rr_radial(rap, la_max, rbp, lb_max, rpw, zetp, eta, mmax, rr(:, :, 1, :))
458 IF (lc_max == 0) RETURN
459 g = 0.5_dp*zetw
460 g_c = 0.5_dp*zetq
461 DO ico = 1, ncoset(max(la_max, lb_max, lc_max))
462 power = indco(:, ico)
463 lowered(:, ico) = [coset(max(0, power(1) - 1), power(2), power(3)), &
464 coset(power(1), max(0, power(2) - 1), power(3)), &
465 coset(power(1), power(2), max(0, power(3) - 1))]
466 END DO
467 DO lc = 1, lc_max
468 DO ic = ncoset(lc - 1) + 1, ncoset(lc)
469 axis = 3
470 IF (indco(3, ic) == 0) axis = 2
471 IF (indco(axis, ic) == 0) axis = 1
472 ic1 = lowered(axis, ic)
473 ic2 = lowered(axis, ic1)
474 DO lb = 0, lb_max
475 DO la = 0, la_max
476 DO m = 0, mmax - la - lb - lc
477 DO ib = ncoset(lb - 1) + 1, ncoset(lb)
478 ib1 = lowered(axis, ib)
479 DO ia = ncoset(la - 1) + 1, ncoset(la)
480 ia1 = lowered(axis, ia)
481 rr(ia, ib, ic, m) = rcw(axis)*rr(ia, ib, ic1, m + 1)
482 IF (indco(axis, ic) > 1) THEN
483 rr(ia, ib, ic, m) = rr(ia, ib, ic, m) + g_c*real(indco(axis, ic) - 1, dp)* &
484 (rr(ia, ib, ic2, m) - eta_c*rr(ia, ib, ic2, m + 1))
485 END IF
486 IF (indco(axis, ia) > 0) THEN
487 rr(ia, ib, ic, m) = rr(ia, ib, ic, m) + &
488 g*real(indco(axis, ia), dp)*rr(ia1, ib, ic1, m + 1)
489 END IF
490 IF (indco(axis, ib) > 0) THEN
491 rr(ia, ib, ic, m) = rr(ia, ib, ic, m) + &
492 g*real(indco(axis, ib), dp)*rr(ia, ib1, ic1, m + 1)
493 END IF
494 END DO
495 END DO
496 END DO
497 END DO
498 END DO
499 END DO
500 END DO
501 END SUBROUTINE os_rr_coulomb3
502
503! **************************************************************************************************
504!> \brief Calculation of the Obara-Saika recurrence relation for 1/r_C
505!> \param rap ...
506!> \param la_max ...
507!> \param rbp ...
508!> \param lb_max ...
509!> \param rcp ...
510!> \param zet ...
511!> \param ldrr1 ...
512!> \param ldrr2 ...
513!> \param rr ...
514!> \date 02.03.2009
515!> \author VW
516!> \version 1.0
517! **************************************************************************************************
518 SUBROUTINE os_rr_coul(rap, la_max, rbp, lb_max, rcp, zet, ldrr1, ldrr2, rr)
519 REAL(dp), DIMENSION(3), INTENT(IN) :: rap
520 INTEGER, INTENT(IN) :: la_max
521 REAL(dp), DIMENSION(3), INTENT(IN) :: rbp
522 INTEGER, INTENT(IN) :: lb_max
523 REAL(dp), DIMENSION(3), INTENT(IN) :: rcp
524 REAL(dp), INTENT(IN) :: zet
525 INTEGER, INTENT(IN) :: ldrr1, ldrr2
526 REAL(dp), DIMENSION(0:ldrr1 - 1, ldrr2, *), &
527 INTENT(INOUT) :: rr
528
529 INTEGER :: mmax
530 REAL(dp) :: t
531
532 mmax = la_max + lb_max
533 t = zet*dot_product(rcp, rcp)
534 CALL fgamma(mmax, t, rr(0:mmax, 1, 1))
535 CALL os_rr_radial_mab(rap, la_max, rbp, lb_max, -rcp, 1.0_dp/zet, 1.0_dp, mmax, &
536 rr(0:mmax, 1:ncoset(la_max), 1:ncoset(lb_max)))
537
538 END SUBROUTINE os_rr_coul
539
540END MODULE ai_os_rr
ai_os_rr::os_rr_ovlp3
Definition ai_os_rr.F:27
ai_os_rr::os_rr_ovlp
Definition ai_os_rr.F:23
ai_os_rr::os_rr_radial
Definition ai_os_rr.F:31
ai_os_rr
Definition ai_os_rr.F:8
ai_os_rr::os_rr_coul
subroutine, public os_rr_coul(rap, la_max, rbp, lb_max, rcp, zet, ldrr1, ldrr2, rr)
Calculation of the Obara-Saika recurrence relation for 1/r_C.
Definition ai_os_rr.F:519
ai_os_rr::os_rr_coulomb3
subroutine, public os_rr_coulomb3(rap, la_max, rbp, lb_max, rcp, lc_max, zetp, zetq, zetw, rr)
Angular three-center Coulomb recurrence from caller-supplied s/s/s seeds.
Definition ai_os_rr.F:443
gamma
Calculation of the incomplete Gamma function F_n(t) for multi-center integrals over Cartesian Gaussia...
Definition gamma.F:15
gamma::fgamma_0
subroutine, public fgamma_0(nmax, t, f)
Calculation of the incomplete Gamma function F(t) for multicenter integrals over Gaussian functions....
Definition gamma.F:154
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