Line data Source code
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 :
8 : ! **************************************************************************************************
9 : !> \brief Contains routines for contraction without dgemms.
10 : !> \notes Contains specific routines for contraction. The compiler flag
11 : !> -D__MAX_CONTR defines the maximum angular momentum up to which
12 : !> specialized code will be used. Default setting is d-functions.
13 : !> Increasing -D__MAX_CONTR produces faster code but might overburden
14 : !> the optimization capabilities of some poor compilers.
15 : !> Templates are available up to l=10. The Fypp parameters lmax_specific
16 : !> and lmax_generic control the generated specialized routines and sparse
17 : !> cases in the generic routine; both default to 4 (g-functions).
18 : !> Fypp preserves compile-time dimensions and explicit sparse sums.
19 : !> \par History
20 : !> 07.2009 created [Manuel Guidon]
21 : !> \author Manuel Guidon
22 : ! **************************************************************************************************
23 :
24 : MODULE hfx_contraction_methods
25 :
26 : !** This defines the default behaviour
27 : #ifndef __MAX_CONTR
28 : #define __MAX_CONTR 2
29 : #endif
30 :
31 : USE kinds, ONLY: dp
32 : #include "../base/base_uses.f90"
33 :
34 : IMPLICIT NONE
35 :
36 : PRIVATE
37 : PUBLIC :: contract
38 :
39 : CONTAINS
40 :
41 : #! Generation limits, independent of the compiler's __MAX_CONTR setting.
42 : #:set lmax_specific = 4
43 : #:set lmax_generic = 4
44 :
45 : #! Shell names and dimensions from the legacy generator, for l=0 through 10.
46 : #:set shells = 'spdfghijklm'
47 : #:set cartesian = [(l + 1)*(l + 2)//2 for l in range(len(shells))]
48 : #:set spherical = [2*l + 1 for l in range(len(shells))]
49 : #:assert 0 <= lmax_specific < len(shells)
50 : #:assert 0 <= lmax_generic < len(shells)
51 : #! Nonzero (Cartesian, spherical) indices, in the original accumulation order.
52 : #:set terms = [ &
53 : [(1, 1)], &
54 : [(1, 3), (2, 1), (3, 2)], &
55 : [(1, 3), (1, 5), (2, 1), (3, 4), (4, 3), (4, 5), (5, 2), (6, 3)], &
56 : [(1, 5), (1, 7), (2, 1), (2, 3), (3, 4), (3, 6), (4, 5), (4, 7), &
57 : (5, 2), (6, 5), (7, 1), (7, 3), (8, 4), (8, 6), (9, 3), (10, 4)], &
58 : [(1, 5), (1, 7), (1, 9), (2, 1), (2, 3), (3, 6), (3, 8), (4, 5), &
59 : (4, 9), (5, 2), (5, 4), (6, 5), (6, 7), (7, 1), (7, 3), (8, 6), &
60 : (8, 8), (9, 3), (10, 6), (11, 5), (11, 7), (11, 9), (12, 2), &
61 : (12, 4), (13, 5), (13, 7), (14, 4), (15, 5)], &
62 : [(1, 7), (1, 9), (1, 11), (2, 1), (2, 3), (2, 5), (3, 6), (3, 8), &
63 : (3, 10), (4, 7), (4, 9), (4, 11), (5, 2), (5, 4), (6, 7), (6, 9), &
64 : (7, 1), (7, 3), (7, 5), (8, 6), (8, 10), (9, 3), (9, 5), (10, 6), &
65 : (10, 8), (11, 7), (11, 9), (11, 11), (12, 2), (12, 4), (13, 7), (13, 9), &
66 : (14, 4), (15, 7), (16, 1), (16, 3), (16, 5), (17, 6), (17, 8), (17, 10), &
67 : (18, 3), (18, 5), (19, 6), (19, 8), (20, 5), (21, 6)], &
68 : [(1, 7), (1, 9), (1, 11), (1, 13), (2, 1), (2, 3), (2, 5), (3, 8), &
69 : (3, 10), (3, 12), (4, 7), (4, 9), (4, 11), (4, 13), (5, 2), (5, 4), &
70 : (5, 6), (6, 7), (6, 9), (6, 11), (7, 1), (7, 5), (8, 8), (8, 10), &
71 : (8, 12), (9, 3), (9, 5), (10, 8), (10, 10), (11, 7), (11, 9), (11, 11), &
72 : (11, 13), (12, 2), (12, 4), (12, 6), (13, 7), (13, 11), (14, 4), (14, 6), &
73 : (15, 7), (15, 9), (16, 1), (16, 3), (16, 5), (17, 8), (17, 10), (17, 12), &
74 : (18, 3), (18, 5), (19, 8), (19, 10), (20, 5), (21, 8), (22, 7), (22, 9), &
75 : (22, 11), (22, 13), (23, 2), (23, 4), (23, 6), (24, 7), (24, 9), (24, 11), &
76 : (25, 4), (25, 6), (26, 7), (26, 9), (27, 6), (28, 7)], &
77 : [(1, 9), (1, 11), (1, 13), (1, 15), (2, 1), (2, 3), (2, 5), (2, 7), &
78 : (3, 8), (3, 10), (3, 12), (3, 14), (4, 9), (4, 11), (4, 13), (4, 15), &
79 : (5, 2), (5, 4), (5, 6), (6, 9), (6, 11), (6, 13), (7, 1), (7, 3), &
80 : (7, 5), (7, 7), (8, 8), (8, 10), (8, 12), (8, 14), (9, 3), (9, 5), &
81 : (9, 7), (10, 8), (10, 10), (10, 12), (11, 9), (11, 11), (11, 13), (11, 15), &
82 : (12, 2), (12, 6), (13, 9), (13, 11), (13, 13), (14, 4), (14, 6), (15, 9), &
83 : (15, 11), (16, 1), (16, 3), (16, 5), (16, 7), (17, 8), (17, 10), (17, 12), &
84 : (17, 14), (18, 3), (18, 5), (18, 7), (19, 8), (19, 12), (20, 5), (20, 7), &
85 : (21, 8), (21, 10), (22, 9), (22, 11), (22, 13), (22, 15), (23, 2), (23, 4), &
86 : (23, 6), (24, 9), (24, 11), (24, 13), (25, 4), (25, 6), (26, 9), (26, 11), &
87 : (27, 6), (28, 9), (29, 1), (29, 3), (29, 5), (29, 7), (30, 8), (30, 10), &
88 : (30, 12), (30, 14), (31, 3), (31, 5), (31, 7), (32, 8), (32, 10), (32, 12), &
89 : (33, 5), (33, 7), (34, 8), (34, 10), (35, 7), (36, 8)], &
90 : [(1, 9), (1, 11), (1, 13), (1, 15), (1, 17), (2, 1), (2, 3), (2, 5), &
91 : (2, 7), (3, 10), (3, 12), (3, 14), (3, 16), (4, 9), (4, 11), (4, 13), &
92 : (4, 15), (4, 17), (5, 2), (5, 4), (5, 6), (5, 8), (6, 9), (6, 11), &
93 : (6, 13), (6, 15), (7, 1), (7, 3), (7, 5), (7, 7), (8, 10), (8, 12), &
94 : (8, 14), (8, 16), (9, 3), (9, 5), (9, 7), (10, 10), (10, 12), (10, 14), &
95 : (11, 9), (11, 13), (11, 17), (12, 2), (12, 4), (12, 6), (12, 8), (13, 9), &
96 : (13, 11), (13, 13), (13, 15), (14, 4), (14, 6), (14, 8), (15, 9), (15, 11), &
97 : (15, 13), (16, 1), (16, 3), (16, 5), (16, 7), (17, 10), (17, 12), (17, 14), &
98 : (17, 16), (18, 3), (18, 7), (19, 10), (19, 12), (19, 14), (20, 5), (20, 7), &
99 : (21, 10), (21, 12), (22, 9), (22, 11), (22, 13), (22, 15), (22, 17), (23, 2), &
100 : (23, 4), (23, 6), (23, 8), (24, 9), (24, 11), (24, 13), (24, 15), (25, 4), &
101 : (25, 6), (25, 8), (26, 9), (26, 13), (27, 6), (27, 8), (28, 9), (28, 11), &
102 : (29, 1), (29, 3), (29, 5), (29, 7), (30, 10), (30, 12), (30, 14), (30, 16), &
103 : (31, 3), (31, 5), (31, 7), (32, 10), (32, 12), (32, 14), (33, 5), (33, 7), &
104 : (34, 10), (34, 12), (35, 7), (36, 10), (37, 9), (37, 11), (37, 13), (37, 15), &
105 : (37, 17), (38, 2), (38, 4), (38, 6), (38, 8), (39, 9), (39, 11), (39, 13), &
106 : (39, 15), (40, 4), (40, 6), (40, 8), (41, 9), (41, 11), (41, 13), (42, 6), &
107 : (42, 8), (43, 9), (43, 11), (44, 8), (45, 9)], &
108 : [(1, 11), (1, 13), (1, 15), (1, 17), (1, 19), (2, 1), (2, 3), (2, 5), &
109 : (2, 7), (2, 9), (3, 10), (3, 12), (3, 14), (3, 16), (3, 18), (4, 11), &
110 : (4, 15), (4, 17), (4, 19), (5, 2), (5, 4), (5, 6), (5, 8), (6, 11), &
111 : (6, 13), (6, 15), (6, 17), (7, 1), (7, 3), (7, 7), (7, 9), (8, 10), &
112 : (8, 12), (8, 14), (8, 16), (8, 18), (9, 3), (9, 5), (9, 7), (9, 9), &
113 : (10, 10), (10, 12), (10, 14), (10, 16), (11, 11), (11, 13), (11, 15), (11, 17), &
114 : (11, 19), (12, 2), (12, 4), (12, 6), (12, 8), (13, 11), (13, 13), (13, 15), &
115 : (13, 17), (14, 4), (14, 6), (14, 8), (15, 11), (15, 13), (15, 15), (16, 1), &
116 : (16, 3), (16, 5), (16, 7), (16, 9), (17, 10), (17, 14), (17, 18), (18, 3), &
117 : (18, 5), (18, 7), (18, 9), (19, 10), (19, 12), (19, 14), (19, 16), (20, 5), &
118 : (20, 7), (20, 9), (21, 10), (21, 12), (21, 14), (22, 11), (22, 13), (22, 17), &
119 : (22, 19), (23, 2), (23, 4), (23, 6), (23, 8), (24, 11), (24, 13), (24, 15), &
120 : (24, 17), (25, 4), (25, 8), (26, 11), (26, 13), (26, 15), (27, 6), (27, 8), &
121 : (28, 11), (28, 13), (29, 1), (29, 3), (29, 5), (29, 9), (30, 10), (30, 12), &
122 : (30, 14), (30, 16), (30, 18), (31, 3), (31, 5), (31, 7), (31, 9), (32, 10), &
123 : (32, 12), (32, 14), (32, 16), (33, 5), (33, 7), (33, 9), (34, 10), (34, 14), &
124 : (35, 7), (35, 9), (36, 10), (36, 12), (37, 11), (37, 13), (37, 15), (37, 17), &
125 : (37, 19), (38, 2), (38, 4), (38, 6), (38, 8), (39, 11), (39, 13), (39, 15), &
126 : (39, 17), (40, 4), (40, 6), (40, 8), (41, 11), (41, 13), (41, 15), (42, 6), &
127 : (42, 8), (43, 11), (43, 13), (44, 8), (45, 11), (46, 1), (46, 3), (46, 5), &
128 : (46, 7), (46, 9), (47, 10), (47, 12), (47, 14), (47, 16), (47, 18), (48, 3), &
129 : (48, 5), (48, 7), (48, 9), (49, 10), (49, 12), (49, 14), (49, 16), (50, 5), &
130 : (50, 7), (50, 9), (51, 10), (51, 12), (51, 14), (52, 7), (52, 9), (53, 10), &
131 : (53, 12), (54, 9), (55, 10)], &
132 : [(1, 11), (1, 13), (1, 15), (1, 17), (1, 19), (1, 21), (2, 1), (2, 3), &
133 : (2, 5), (2, 7), (2, 9), (3, 12), (3, 14), (3, 16), (3, 18), (3, 20), &
134 : (4, 11), (4, 13), (4, 15), (4, 17), (4, 19), (4, 21), (5, 2), (5, 4), &
135 : (5, 6), (5, 8), (5, 10), (6, 11), (6, 13), (6, 15), (6, 17), (6, 19), &
136 : (7, 1), (7, 3), (7, 5), (7, 7), (7, 9), (8, 12), (8, 16), (8, 18), &
137 : (8, 20), (9, 3), (9, 5), (9, 7), (9, 9), (10, 12), (10, 14), (10, 16), &
138 : (10, 18), (11, 11), (11, 13), (11, 15), (11, 17), (11, 19), (11, 21), (12, 2), &
139 : (12, 4), (12, 8), (12, 10), (13, 11), (13, 13), (13, 15), (13, 17), (13, 19), &
140 : (14, 4), (14, 6), (14, 8), (14, 10), (15, 11), (15, 13), (15, 15), (15, 17), &
141 : (16, 1), (16, 5), (16, 9), (17, 12), (17, 14), (17, 16), (17, 18), (17, 20), &
142 : (18, 3), (18, 5), (18, 7), (18, 9), (19, 12), (19, 14), (19, 16), (19, 18), &
143 : (20, 5), (20, 7), (20, 9), (21, 12), (21, 14), (21, 16), (22, 11), (22, 13), &
144 : (22, 15), (22, 17), (22, 19), (22, 21), (23, 2), (23, 4), (23, 6), (23, 8), &
145 : (23, 10), (24, 11), (24, 15), (24, 19), (25, 4), (25, 6), (25, 8), (25, 10), &
146 : (26, 11), (26, 13), (26, 15), (26, 17), (27, 6), (27, 8), (27, 10), (28, 11), &
147 : (28, 13), (28, 15), (29, 1), (29, 3), (29, 5), (29, 7), (29, 9), (30, 12), &
148 : (30, 14), (30, 18), (30, 20), (31, 3), (31, 5), (31, 7), (31, 9), (32, 12), &
149 : (32, 14), (32, 16), (32, 18), (33, 5), (33, 9), (34, 12), (34, 14), (34, 16), &
150 : (35, 7), (35, 9), (36, 12), (36, 14), (37, 11), (37, 13), (37, 15), (37, 17), &
151 : (37, 19), (37, 21), (38, 2), (38, 4), (38, 6), (38, 10), (39, 11), (39, 13), &
152 : (39, 15), (39, 17), (39, 19), (40, 4), (40, 6), (40, 8), (40, 10), (41, 11), &
153 : (41, 13), (41, 15), (41, 17), (42, 6), (42, 8), (42, 10), (43, 11), (43, 15), &
154 : (44, 8), (44, 10), (45, 11), (45, 13), (46, 1), (46, 3), (46, 5), (46, 7), &
155 : (46, 9), (47, 12), (47, 14), (47, 16), (47, 18), (47, 20), (48, 3), (48, 5), &
156 : (48, 7), (48, 9), (49, 12), (49, 14), (49, 16), (49, 18), (50, 5), (50, 7), &
157 : (50, 9), (51, 12), (51, 14), (51, 16), (52, 7), (52, 9), (53, 12), (53, 14), &
158 : (54, 9), (55, 12), (56, 11), (56, 13), (56, 15), (56, 17), (56, 19), (56, 21), &
159 : (57, 2), (57, 4), (57, 6), (57, 8), (57, 10), (58, 11), (58, 13), (58, 15), &
160 : (58, 17), (58, 19), (59, 4), (59, 6), (59, 8), (59, 10), (60, 11), (60, 13), &
161 : (60, 15), (60, 17), (61, 6), (61, 8), (61, 10), (62, 11), (62, 13), (62, 15), &
162 : (63, 8), (63, 10), (64, 11), (64, 13), (65, 10), (66, 11)]]
163 :
164 : #:def call_generic()
165 : CALL contract_generic(ncoa, ncob, ncoc, ncod, nsoa, nsob, nsoc, nsod, &
166 : work, nl_a, nl_b, nl_c, nl_d, &
167 : sphi_a, sphi_b, sphi_c, sphi_d, primitives, buffer1, buffer2)
168 : #:enddef
169 :
170 : #! The first three stages transpose into a work buffer; the last accumulates
171 : #! into the four-dimensional output. All sparse terms are expanded by Fypp.
172 : #:def transform(angular, stage, ns)
173 : #:set axis = 'abcd'[stage]
174 : #:set source = ['work', 'buffer1', 'buffer2', 'buffer1'][stage]
175 : #:if stage < 3
176 : #:set target = ['buffer1', 'buffer2', 'buffer1'][stage]
177 : DO i = 1, imax
178 : #:for k, j in terms[angular]
179 : ${target}$ (i + imax*(${j}$-1)) = ${target}$ (i + imax*(${j}$-1)) &
180 : + ${source}$ (${k}$+(i - 1)*kmax)*sphi_${axis}$ (${k}$, ${j}$+s_offset_${axis}$1)
181 : #:endfor
182 : END DO
183 : #:else
184 : DO i1 = 1, ${ns[2]}$
185 : DO i2 = 1, ${ns[1]}$
186 : DO i3 = 1, ${ns[0]}$
187 : i = i + 1
188 : #:for k, j in terms[angular]
189 : primitives(s_offset_a1 + i3, s_offset_b1 + i2, s_offset_c1 + i1, s_offset_d1 + ${j}$) = &
190 : primitives(s_offset_a1 + i3, s_offset_b1 + i2, s_offset_c1 + i1, s_offset_d1 + ${j}$) &
191 : + buffer1(${k}$+(i - 1)*kmax)*sphi_d(${k}$, ${j}$+s_offset_d1)
192 : #:endfor
193 : END DO
194 : END DO
195 : END DO
196 : #:endif
197 : #:enddef
198 :
199 : #! Share the loop structure between all specialized kernels and the fallback.
200 : #! Runtime dimensions and dense fallback loops are emitted only for generic.
201 : #:def contraction(angular=None)
202 : #:set generic = angular is None
203 : #:if generic
204 : #:set name = 'generic'
205 : #:set nc = ['nco' + axis for axis in 'abcd']
206 : #:set ns = ['nso' + axis for axis in 'abcd']
207 : #:else
208 : #:set name = ''.join(shells[l] for l in angular)
209 : #:set nc = [str(cartesian[l]) for l in angular]
210 : #:set ns = [str(spherical[l]) for l in angular]
211 : #:set max_contr = max(angular)
212 : #if __MAX_CONTR >= ${max_contr}$
213 : #:endif
214 : ! **************************************************************************************************
215 : !> \brief Contract Cartesian integrals into spherical basis functions.
216 : #:if generic
217 : #:for axis in 'abcd'
218 : !> \param nco${axis}$ Number of Cartesian functions on center ${axis}$
219 : #:endfor
220 : #:for axis in 'abcd'
221 : !> \param nso${axis}$ Number of spherical functions on center ${axis}$
222 : #:endfor
223 : #:endif
224 : !> \param work Cartesian integrals
225 : #:for axis in 'abcd'
226 : !> \param nl_${axis}$ Number of contractions on center ${axis}$
227 : #:endfor
228 : #:for axis in 'abcd'
229 : !> \param sphi_${axis}$ Contraction and spherical transformation coefficients
230 : #:endfor
231 : !> \param primitives Accumulated spherical integrals
232 : !> \param buffer1 Work buffer
233 : !> \param buffer2 Work buffer
234 : ! **************************************************************************************************
235 : #:if generic
236 : SUBROUTINE contract_generic(ncoa, ncob, ncoc, ncod, nsoa, nsob, nsoc, nsod, &
237 : work, nl_a, nl_b, nl_c, nl_d, &
238 : sphi_a, sphi_b, sphi_c, sphi_d, primitives, buffer1, buffer2)
239 :
240 : INTEGER, INTENT(IN) :: ncoa, ncob, ncoc, ncod, nsoa, nsob, nsoc, nsod
241 : #:else
242 : SUBROUTINE contract_${name}$ (work, nl_a, nl_b, nl_c, nl_d, &
243 : sphi_a, sphi_b, sphi_c, sphi_d, primitives, buffer1, buffer2)
244 :
245 : #:endif
246 : REAL(dp), DIMENSION(${ '*'.join(nc) }$), INTENT(IN) :: work
247 : INTEGER :: nl_a, nl_b, nl_c, nl_d
248 : #:for axis, ncart, nsph in zip('abcd', nc, ns)
249 : REAL(dp), DIMENSION(${ncart}$, ${nsph}$*nl_${axis}$), INTENT(IN) :: sphi_${axis}$
250 : #:endfor
251 : REAL(dp), DIMENSION(${ns[0]}$*nl_a, ${ns[1]}$*nl_b, ${ns[2]}$*nl_c, ${ns[3]}$*nl_d) :: primitives
252 : REAL(dp), DIMENSION(${ '*'.join(nc) }$) :: buffer1, buffer2
253 :
254 : INTEGER :: i, i1, i2, i3, ia, ib, ic, id, imax, &
255 : #:if generic
256 : j, jmax, k, &
257 : #:endif
258 : kmax, s_offset_a1, s_offset_b1, s_offset_c1, s_offset_d1
259 :
260 : s_offset_a1 = 0
261 : DO ia = 1, nl_a
262 : s_offset_b1 = 0
263 : DO ib = 1, nl_b
264 : s_offset_c1 = 0
265 : DO ic = 1, nl_c
266 : s_offset_d1 = 0
267 : DO id = 1, nl_d
268 : #:for stage, axis in enumerate('abcd')
269 : #:set source = ['work', 'buffer1', 'buffer2', 'buffer1'][stage]
270 : #:if stage < 3
271 : #:set target = ['buffer1', 'buffer2', 'buffer1'][stage]
272 : ${target}$ = 0.0_dp
273 : #:endif
274 : imax = ${ '*'.join(ns[:stage] + nc[stage + 1:]) }$
275 : #:if generic
276 : jmax = ${ns[stage]}$
277 : #:endif
278 : kmax = ${nc[stage]}$
279 : #:if stage == 3
280 : i = 0
281 : #:endif
282 : #:if generic
283 : SELECT CASE (kmax)
284 : #:for l, ncart in enumerate(cartesian[:lmax_generic + 1])
285 : CASE (${ncart}$)
286 : $:transform(l, stage, ns)
287 : #:endfor
288 : CASE DEFAULT
289 : DO j = 1, jmax
290 : #:if stage < 3
291 : DO i = 1, imax
292 : DO k = 1, kmax
293 : ${target}$ (i + imax*(j - 1)) = ${target}$ (i + imax*(j - 1)) &
294 : + ${source}$ (k + (i - 1)*kmax)*sphi_${axis}$ (k, j + s_offset_${axis}$1)
295 : END DO
296 : END DO
297 : #:else
298 : i = 0
299 : DO i1 = 1, nsoc
300 : DO i2 = 1, nsob
301 : DO i3 = 1, nsoa
302 : i = i + 1
303 : DO k = 1, kmax
304 : primitives(s_offset_a1 + i3, s_offset_b1 + i2, s_offset_c1 + i1, s_offset_d1 + j) = &
305 : primitives(s_offset_a1 + i3, s_offset_b1 + i2, s_offset_c1 + i1, s_offset_d1 + j) &
306 : + buffer1(k + (i - 1)*kmax)*sphi_d(k, j + s_offset_d1)
307 : END DO
308 : END DO
309 : END DO
310 : END DO
311 : #:endif
312 : END DO
313 : END SELECT
314 : #:else
315 : $:transform(angular[stage], stage, ns)
316 : #:endif
317 : #:endfor
318 : s_offset_d1 = s_offset_d1 + ${ns[3]}$
319 : END DO
320 : s_offset_c1 = s_offset_c1 + ${ns[2]}$
321 : END DO
322 : s_offset_b1 = s_offset_b1 + ${ns[1]}$
323 : END DO
324 : s_offset_a1 = s_offset_a1 + ${ns[0]}$
325 : END DO
326 : END SUBROUTINE contract_${name}$
327 : #:if not generic
328 : #endif
329 : #:endif
330 : #:enddef
331 :
332 : ! **************************************************************************************************
333 : !> \brief ...
334 : !> \param ncoa ...
335 : !> \param ncob ...
336 : !> \param ncoc ...
337 : !> \param ncod ...
338 : !> \param nsoa ...
339 : !> \param nsob ...
340 : !> \param nsoc ...
341 : !> \param nsod ...
342 : !> \param n_a ...
343 : !> \param n_b ...
344 : !> \param n_c ...
345 : !> \param n_d ...
346 : !> \param nl_a ...
347 : !> \param nl_b ...
348 : !> \param nl_c ...
349 : !> \param nl_d ...
350 : !> \param work ...
351 : !> \param sphi_a ...
352 : !> \param sphi_b ...
353 : !> \param sphi_c ...
354 : !> \param sphi_d ...
355 : !> \param primitives ...
356 : !> \param buffer1 ...
357 : !> \param buffer2 ...
358 : ! **************************************************************************************************
359 526124097 : SUBROUTINE contract(ncoa, ncob, ncoc, ncod, nsoa, nsob, nsoc, nsod, &
360 526124097 : n_a, n_b, n_c, n_d, nl_a, nl_b, nl_c, nl_d, work, &
361 526124097 : sphi_a, sphi_b, sphi_c, sphi_d, &
362 526124097 : primitives, &
363 526124097 : buffer1, buffer2)
364 :
365 : INTEGER, INTENT(IN) :: ncoa, ncob, ncoc, ncod, nsoa, nsob, nsoc, nsod, &
366 : n_a, n_b, n_c, n_d, nl_a, nl_b, nl_c, nl_d
367 : REAL(dp), DIMENSION(ncoa*ncob*ncoc*ncod), INTENT(IN) :: work
368 : REAL(dp), DIMENSION(ncoa, nsoa*nl_a), INTENT(IN) :: sphi_a
369 : REAL(dp), DIMENSION(ncob, nsob*nl_b), INTENT(IN) :: sphi_b
370 : REAL(dp), DIMENSION(ncoc, nsoc*nl_c), INTENT(IN) :: sphi_c
371 : REAL(dp), DIMENSION(ncod, nsod*nl_d), INTENT(IN) :: sphi_d
372 :
373 : REAL(dp), DIMENSION(nsoa*nl_a, nsob*nl_b, nsoc*nl_c, nsod*nl_d) :: primitives
374 : REAL(dp), DIMENSION(ncoa*ncob*ncoc*ncod) :: buffer1, buffer2
375 :
376 : #if !defined(__LIBINT)
377 : MARK_USED(ncoa)
378 : MARK_USED(ncob)
379 : MARK_USED(ncoc)
380 : MARK_USED(ncod)
381 : MARK_USED(nsoa)
382 : MARK_USED(nsob)
383 : MARK_USED(nsoc)
384 : MARK_USED(nsod)
385 : MARK_USED(n_a)
386 : MARK_USED(n_b)
387 : MARK_USED(n_c)
388 : MARK_USED(n_d)
389 : MARK_USED(nl_a)
390 : MARK_USED(nl_b)
391 : MARK_USED(nl_c)
392 : MARK_USED(nl_d)
393 : MARK_USED(sphi_a)
394 : MARK_USED(sphi_b)
395 : MARK_USED(sphi_c)
396 : MARK_USED(sphi_d)
397 : MARK_USED(work)
398 : MARK_USED(primitives)
399 : MARK_USED(buffer1)
400 : MARK_USED(buffer2)
401 : CPABORT("libint not compiled in")
402 : #else
403 875634164 : SELECT CASE (n_a)
404 : #:for a in range(lmax_specific + 1)
405 : CASE (${a}$)
406 949479850 : SELECT CASE (n_b)
407 : #:for b in range(lmax_specific + 1)
408 : CASE (${b}$)
409 823668048 : SELECT CASE (n_c)
410 : #:for c in range(lmax_specific + 1)
411 : CASE (${c}$)
412 882746491 : SELECT CASE (n_d)
413 : #:for d in range(lmax_specific + 1)
414 : CASE (${d}$)
415 : #:set max_contr = max(a, b, c, d)
416 : #if __MAX_CONTR >= ${max_contr}$
417 : CALL contract_${''.join(shells[l] for l in (a, b, c, d))}$ (work, &
418 : nl_a, nl_b, nl_c, nl_d, &
419 525343160 : sphi_a, sphi_b, sphi_c, sphi_d, primitives, buffer1, buffer2)
420 : #else
421 780937 : $:call_generic()
422 : #endif
423 : #:endfor
424 : CASE DEFAULT
425 526124097 : $:call_generic()
426 : END SELECT
427 : #:endfor
428 : CASE DEFAULT
429 526124097 : $:call_generic()
430 : END SELECT
431 : #:endfor
432 : CASE DEFAULT
433 526124097 : $:call_generic()
434 : END SELECT
435 : #:endfor
436 : CASE DEFAULT
437 526124097 : $:call_generic()
438 : END SELECT
439 : #endif
440 526124097 : END SUBROUTINE contract
441 :
442 : #if defined(__LIBINT)
443 : #:for a in range(lmax_specific + 1)
444 : #:for b in range(lmax_specific + 1)
445 : #:for c in range(lmax_specific + 1)
446 : #:for d in range(lmax_specific + 1)
447 64061089671 : $:contraction((a, b, c, d))
448 : #:endfor
449 : #:endfor
450 : #:endfor
451 : #:endfor
452 :
453 711974386 : $:contraction()
454 : #endif
455 : END MODULE hfx_contraction_methods
|
