1 #include "simint/boys/boys.h"
2 #include "simint/ostei/gen/ostei_generated.h"
3 #include "simint/vectorization/vectorization.h"
4 #include <math.h>
5 #include <string.h>
6
7
ostei_f_p_i_h(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__f_p_i_h)8 int ostei_f_p_i_h(struct simint_multi_shellpair const P,
9 struct simint_multi_shellpair const Q,
10 double screen_tol,
11 double * const restrict work,
12 double * const restrict INT__f_p_i_h)
13 {
14
15 SIMINT_ASSUME_ALIGN_DBL(work);
16 SIMINT_ASSUME_ALIGN_DBL(INT__f_p_i_h);
17 int ab, cd, abcd;
18 int istart, jstart;
19 int iprimcd, nprim_icd, icd;
20 const int check_screen = (screen_tol > 0.0);
21 int i, j;
22 int n;
23 int not_screened;
24 int real_abcd;
25 int iket;
26 int ibra;
27
28 // partition workspace
29 double * const INT__f_s_i_s = work + (SIMINT_NSHELL_SIMD * 0);
30 double * const INT__f_s_k_s = work + (SIMINT_NSHELL_SIMD * 280);
31 double * const INT__f_s_l_s = work + (SIMINT_NSHELL_SIMD * 640);
32 double * const INT__f_s_m_s = work + (SIMINT_NSHELL_SIMD * 1090);
33 double * const INT__f_s_n_s = work + (SIMINT_NSHELL_SIMD * 1640);
34 double * const INT__f_s_o_s = work + (SIMINT_NSHELL_SIMD * 2300);
35 double * const INT__g_s_i_s = work + (SIMINT_NSHELL_SIMD * 3080);
36 double * const INT__g_s_k_s = work + (SIMINT_NSHELL_SIMD * 3500);
37 double * const INT__g_s_l_s = work + (SIMINT_NSHELL_SIMD * 4040);
38 double * const INT__g_s_m_s = work + (SIMINT_NSHELL_SIMD * 4715);
39 double * const INT__g_s_n_s = work + (SIMINT_NSHELL_SIMD * 5540);
40 double * const INT__g_s_o_s = work + (SIMINT_NSHELL_SIMD * 6530);
41 SIMINT_DBLTYPE * const primwork = (SIMINT_DBLTYPE *)(work + SIMINT_NSHELL_SIMD*7700);
42 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_s = primwork + 0;
43 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_p_s = primwork + 16;
44 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_d_s = primwork + 61;
45 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_f_s = primwork + 145;
46 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_g_s = primwork + 275;
47 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_h_s = primwork + 455;
48 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_i_s = primwork + 686;
49 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_k_s = primwork + 966;
50 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_l_s = primwork + 1290;
51 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_m_s = primwork + 1650;
52 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_n_s = primwork + 2035;
53 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_o_s = primwork + 2431;
54 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_f_s = primwork + 2821;
55 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_g_s = primwork + 2941;
56 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_h_s = primwork + 3121;
57 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_i_s = primwork + 3373;
58 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_k_s = primwork + 3709;
59 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_l_s = primwork + 4141;
60 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_m_s = primwork + 4681;
61 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_n_s = primwork + 5341;
62 SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_o_s = primwork + 6133;
63 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_g_s = primwork + 7069;
64 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_h_s = primwork + 7339;
65 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_i_s = primwork + 7717;
66 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_k_s = primwork + 8221;
67 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_l_s = primwork + 8869;
68 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_m_s = primwork + 9679;
69 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_n_s = primwork + 10669;
70 SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_o_s = primwork + 11857;
71 SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_h_s = primwork + 13261;
72 SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_i_s = primwork + 13681;
73 SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_k_s = primwork + 14241;
74 SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_l_s = primwork + 14961;
75 SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_m_s = primwork + 15861;
76 SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_n_s = primwork + 16961;
77 SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_o_s = primwork + 18281;
78 SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_i_s = primwork + 19841;
79 SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_k_s = primwork + 20261;
80 SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_l_s = primwork + 20801;
81 SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_m_s = primwork + 21476;
82 SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_n_s = primwork + 22301;
83 SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_o_s = primwork + 23291;
84 double * const hrrwork = (double *)(primwork + 24461);
85 double * const HRR_INT__f_p_i_s = hrrwork + 0;
86 double * const HRR_INT__f_p_i_p = hrrwork + 840;
87 double * const HRR_INT__f_p_i_d = hrrwork + 3360;
88 double * const HRR_INT__f_p_i_f = hrrwork + 8400;
89 double * const HRR_INT__f_p_i_g = hrrwork + 16800;
90 double * const HRR_INT__f_p_k_s = hrrwork + 29400;
91 double * const HRR_INT__f_p_k_p = hrrwork + 30480;
92 double * const HRR_INT__f_p_k_d = hrrwork + 33720;
93 double * const HRR_INT__f_p_k_f = hrrwork + 40200;
94 double * const HRR_INT__f_p_k_g = hrrwork + 51000;
95 double * const HRR_INT__f_p_l_s = hrrwork + 67200;
96 double * const HRR_INT__f_p_l_p = hrrwork + 68550;
97 double * const HRR_INT__f_p_l_d = hrrwork + 72600;
98 double * const HRR_INT__f_p_l_f = hrrwork + 80700;
99 double * const HRR_INT__f_p_m_s = hrrwork + 94200;
100 double * const HRR_INT__f_p_m_p = hrrwork + 95850;
101 double * const HRR_INT__f_p_m_d = hrrwork + 100800;
102 double * const HRR_INT__f_p_n_s = hrrwork + 110700;
103 double * const HRR_INT__f_p_n_p = hrrwork + 112680;
104 double * const HRR_INT__f_p_o_s = hrrwork + 118620;
105
106
107 // Create constants
108 const SIMINT_DBLTYPE const_1 = SIMINT_DBLSET1(1);
109 const SIMINT_DBLTYPE const_10 = SIMINT_DBLSET1(10);
110 const SIMINT_DBLTYPE const_11 = SIMINT_DBLSET1(11);
111 const SIMINT_DBLTYPE const_2 = SIMINT_DBLSET1(2);
112 const SIMINT_DBLTYPE const_3 = SIMINT_DBLSET1(3);
113 const SIMINT_DBLTYPE const_4 = SIMINT_DBLSET1(4);
114 const SIMINT_DBLTYPE const_5 = SIMINT_DBLSET1(5);
115 const SIMINT_DBLTYPE const_6 = SIMINT_DBLSET1(6);
116 const SIMINT_DBLTYPE const_7 = SIMINT_DBLSET1(7);
117 const SIMINT_DBLTYPE const_8 = SIMINT_DBLSET1(8);
118 const SIMINT_DBLTYPE const_9 = SIMINT_DBLSET1(9);
119 const SIMINT_DBLTYPE one_half = SIMINT_DBLSET1(0.5);
120
121
122 ////////////////////////////////////////
123 // Loop over shells and primitives
124 ////////////////////////////////////////
125
126 real_abcd = 0;
127 istart = 0;
128 for(ab = 0; ab < P.nshell12_clip; ++ab)
129 {
130 const int iend = istart + P.nprim12[ab];
131
132 cd = 0;
133 jstart = 0;
134
135 for(cd = 0; cd < Q.nshell12_clip; cd += SIMINT_NSHELL_SIMD)
136 {
137 const int nshellbatch = ((cd + SIMINT_NSHELL_SIMD) > Q.nshell12_clip) ? Q.nshell12_clip - cd : SIMINT_NSHELL_SIMD;
138 int jend = jstart;
139 for(i = 0; i < nshellbatch; i++)
140 jend += Q.nprim12[cd+i];
141
142 // Clear the beginning of the workspace (where we are accumulating integrals)
143 memset(work, 0, SIMINT_NSHELL_SIMD * 7700 * sizeof(double));
144 abcd = 0;
145
146
147 for(i = istart; i < iend; ++i)
148 {
149 SIMINT_DBLTYPE bra_screen_max; // only used if check_screen
150
151 if(check_screen)
152 {
153 // Skip this whole thing if always insignificant
154 if((P.screen[i] * Q.screen_max) < screen_tol)
155 continue;
156 bra_screen_max = SIMINT_DBLSET1(P.screen[i]);
157 }
158
159 icd = 0;
160 iprimcd = 0;
161 nprim_icd = Q.nprim12[cd];
162 double * restrict PRIM_PTR_INT__f_s_i_s = INT__f_s_i_s + abcd * 280;
163 double * restrict PRIM_PTR_INT__f_s_k_s = INT__f_s_k_s + abcd * 360;
164 double * restrict PRIM_PTR_INT__f_s_l_s = INT__f_s_l_s + abcd * 450;
165 double * restrict PRIM_PTR_INT__f_s_m_s = INT__f_s_m_s + abcd * 550;
166 double * restrict PRIM_PTR_INT__f_s_n_s = INT__f_s_n_s + abcd * 660;
167 double * restrict PRIM_PTR_INT__f_s_o_s = INT__f_s_o_s + abcd * 780;
168 double * restrict PRIM_PTR_INT__g_s_i_s = INT__g_s_i_s + abcd * 420;
169 double * restrict PRIM_PTR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
170 double * restrict PRIM_PTR_INT__g_s_l_s = INT__g_s_l_s + abcd * 675;
171 double * restrict PRIM_PTR_INT__g_s_m_s = INT__g_s_m_s + abcd * 825;
172 double * restrict PRIM_PTR_INT__g_s_n_s = INT__g_s_n_s + abcd * 990;
173 double * restrict PRIM_PTR_INT__g_s_o_s = INT__g_s_o_s + abcd * 1170;
174
175
176
177 // Load these one per loop over i
178 const SIMINT_DBLTYPE P_alpha = SIMINT_DBLSET1(P.alpha[i]);
179 const SIMINT_DBLTYPE P_prefac = SIMINT_DBLSET1(P.prefac[i]);
180 const SIMINT_DBLTYPE Pxyz[3] = { SIMINT_DBLSET1(P.x[i]), SIMINT_DBLSET1(P.y[i]), SIMINT_DBLSET1(P.z[i]) };
181
182 const SIMINT_DBLTYPE P_PA[3] = { SIMINT_DBLSET1(P.PA_x[i]), SIMINT_DBLSET1(P.PA_y[i]), SIMINT_DBLSET1(P.PA_z[i]) };
183
184 for(j = jstart; j < jend; j += SIMINT_SIMD_LEN)
185 {
186 // calculate the shell offsets
187 // these are the offset from the shell pointed to by cd
188 // for each element
189 int shelloffsets[SIMINT_SIMD_LEN] = {0};
190 int lastoffset = 0;
191 const int nlane = ( ((j + SIMINT_SIMD_LEN) < jend) ? SIMINT_SIMD_LEN : (jend - j));
192
193 if((iprimcd + SIMINT_SIMD_LEN) >= nprim_icd)
194 {
195 // Handle if the first element of the vector is a new shell
196 if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
197 {
198 nprim_icd += Q.nprim12[cd + (++icd)];
199 PRIM_PTR_INT__f_s_i_s += 280;
200 PRIM_PTR_INT__f_s_k_s += 360;
201 PRIM_PTR_INT__f_s_l_s += 450;
202 PRIM_PTR_INT__f_s_m_s += 550;
203 PRIM_PTR_INT__f_s_n_s += 660;
204 PRIM_PTR_INT__f_s_o_s += 780;
205 PRIM_PTR_INT__g_s_i_s += 420;
206 PRIM_PTR_INT__g_s_k_s += 540;
207 PRIM_PTR_INT__g_s_l_s += 675;
208 PRIM_PTR_INT__g_s_m_s += 825;
209 PRIM_PTR_INT__g_s_n_s += 990;
210 PRIM_PTR_INT__g_s_o_s += 1170;
211 }
212 iprimcd++;
213 for(n = 1; n < SIMINT_SIMD_LEN; ++n)
214 {
215 if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
216 {
217 shelloffsets[n] = shelloffsets[n-1] + 1;
218 lastoffset++;
219 nprim_icd += Q.nprim12[cd + (++icd)];
220 }
221 else
222 shelloffsets[n] = shelloffsets[n-1];
223 iprimcd++;
224 }
225 }
226 else
227 iprimcd += SIMINT_SIMD_LEN;
228
229 // Do we have to compute this vector (or has it been screened out)?
230 // (not_screened != 0 means we have to do this vector)
231 if(check_screen)
232 {
233 const double vmax = vector_max(SIMINT_MUL(bra_screen_max, SIMINT_DBLLOAD(Q.screen, j)));
234 if(vmax < screen_tol)
235 {
236 PRIM_PTR_INT__f_s_i_s += lastoffset*280;
237 PRIM_PTR_INT__f_s_k_s += lastoffset*360;
238 PRIM_PTR_INT__f_s_l_s += lastoffset*450;
239 PRIM_PTR_INT__f_s_m_s += lastoffset*550;
240 PRIM_PTR_INT__f_s_n_s += lastoffset*660;
241 PRIM_PTR_INT__f_s_o_s += lastoffset*780;
242 PRIM_PTR_INT__g_s_i_s += lastoffset*420;
243 PRIM_PTR_INT__g_s_k_s += lastoffset*540;
244 PRIM_PTR_INT__g_s_l_s += lastoffset*675;
245 PRIM_PTR_INT__g_s_m_s += lastoffset*825;
246 PRIM_PTR_INT__g_s_n_s += lastoffset*990;
247 PRIM_PTR_INT__g_s_o_s += lastoffset*1170;
248 continue;
249 }
250 }
251
252 const SIMINT_DBLTYPE Q_alpha = SIMINT_DBLLOAD(Q.alpha, j);
253 const SIMINT_DBLTYPE PQalpha_mul = SIMINT_MUL(P_alpha, Q_alpha);
254 const SIMINT_DBLTYPE PQalpha_sum = SIMINT_ADD(P_alpha, Q_alpha);
255 const SIMINT_DBLTYPE one_over_PQalpha_sum = SIMINT_DIV(const_1, PQalpha_sum);
256
257
258 /* construct R2 = (Px - Qx)**2 + (Py - Qy)**2 + (Pz -Qz)**2 */
259 SIMINT_DBLTYPE PQ[3];
260 PQ[0] = SIMINT_SUB(Pxyz[0], SIMINT_DBLLOAD(Q.x, j));
261 PQ[1] = SIMINT_SUB(Pxyz[1], SIMINT_DBLLOAD(Q.y, j));
262 PQ[2] = SIMINT_SUB(Pxyz[2], SIMINT_DBLLOAD(Q.z, j));
263 SIMINT_DBLTYPE R2 = SIMINT_MUL(PQ[0], PQ[0]);
264 R2 = SIMINT_FMADD(PQ[1], PQ[1], R2);
265 R2 = SIMINT_FMADD(PQ[2], PQ[2], R2);
266
267 const SIMINT_DBLTYPE alpha = SIMINT_MUL(PQalpha_mul, one_over_PQalpha_sum); // alpha from MEST
268 const SIMINT_DBLTYPE one_over_p = SIMINT_DIV(const_1, P_alpha);
269 const SIMINT_DBLTYPE one_over_q = SIMINT_DIV(const_1, Q_alpha);
270 const SIMINT_DBLTYPE one_over_2p = SIMINT_MUL(one_half, one_over_p);
271 const SIMINT_DBLTYPE one_over_2q = SIMINT_MUL(one_half, one_over_q);
272 const SIMINT_DBLTYPE one_over_2pq = SIMINT_MUL(one_half, one_over_PQalpha_sum);
273 const SIMINT_DBLTYPE Q_PA[3] = { SIMINT_DBLLOAD(Q.PA_x, j), SIMINT_DBLLOAD(Q.PA_y, j), SIMINT_DBLLOAD(Q.PA_z, j) };
274
275 // NOTE: Minus sign!
276 const SIMINT_DBLTYPE a_over_p = SIMINT_MUL(SIMINT_NEG(alpha), one_over_p);
277 SIMINT_DBLTYPE aop_PQ[3];
278 aop_PQ[0] = SIMINT_MUL(a_over_p, PQ[0]);
279 aop_PQ[1] = SIMINT_MUL(a_over_p, PQ[1]);
280 aop_PQ[2] = SIMINT_MUL(a_over_p, PQ[2]);
281
282 SIMINT_DBLTYPE a_over_q = SIMINT_MUL(alpha, one_over_q);
283 SIMINT_DBLTYPE aoq_PQ[3];
284 aoq_PQ[0] = SIMINT_MUL(a_over_q, PQ[0]);
285 aoq_PQ[1] = SIMINT_MUL(a_over_q, PQ[1]);
286 aoq_PQ[2] = SIMINT_MUL(a_over_q, PQ[2]);
287 // Put a minus sign here so we don't have to in RR routines
288 a_over_q = SIMINT_NEG(a_over_q);
289
290
291 //////////////////////////////////////////////
292 // Fjt function section
293 // Maximum v value: 15
294 //////////////////////////////////////////////
295 // The parameter to the Fjt function
296 const SIMINT_DBLTYPE F_x = SIMINT_MUL(R2, alpha);
297
298
299 const SIMINT_DBLTYPE Q_prefac = mask_load(nlane, Q.prefac + j);
300
301
302 boys_F_split(PRIM_INT__s_s_s_s, F_x, 15);
303 SIMINT_DBLTYPE prefac = SIMINT_SQRT(one_over_PQalpha_sum);
304 prefac = SIMINT_MUL(SIMINT_MUL(P_prefac, Q_prefac), prefac);
305 for(n = 0; n <= 15; n++)
306 PRIM_INT__s_s_s_s[n] = SIMINT_MUL(PRIM_INT__s_s_s_s[n], prefac);
307
308 //////////////////////////////////////////////
309 // Primitive integrals: Vertical recurrance
310 //////////////////////////////////////////////
311
312 const SIMINT_DBLTYPE vrr_const_1_over_2p = one_over_2p;
313 const SIMINT_DBLTYPE vrr_const_2_over_2p = SIMINT_MUL(const_2, one_over_2p);
314 const SIMINT_DBLTYPE vrr_const_3_over_2p = SIMINT_MUL(const_3, one_over_2p);
315 const SIMINT_DBLTYPE vrr_const_1_over_2q = one_over_2q;
316 const SIMINT_DBLTYPE vrr_const_2_over_2q = SIMINT_MUL(const_2, one_over_2q);
317 const SIMINT_DBLTYPE vrr_const_3_over_2q = SIMINT_MUL(const_3, one_over_2q);
318 const SIMINT_DBLTYPE vrr_const_4_over_2q = SIMINT_MUL(const_4, one_over_2q);
319 const SIMINT_DBLTYPE vrr_const_5_over_2q = SIMINT_MUL(const_5, one_over_2q);
320 const SIMINT_DBLTYPE vrr_const_6_over_2q = SIMINT_MUL(const_6, one_over_2q);
321 const SIMINT_DBLTYPE vrr_const_7_over_2q = SIMINT_MUL(const_7, one_over_2q);
322 const SIMINT_DBLTYPE vrr_const_8_over_2q = SIMINT_MUL(const_8, one_over_2q);
323 const SIMINT_DBLTYPE vrr_const_9_over_2q = SIMINT_MUL(const_9, one_over_2q);
324 const SIMINT_DBLTYPE vrr_const_10_over_2q = SIMINT_MUL(const_10, one_over_2q);
325 const SIMINT_DBLTYPE vrr_const_1_over_2pq = one_over_2pq;
326 const SIMINT_DBLTYPE vrr_const_2_over_2pq = SIMINT_MUL(const_2, one_over_2pq);
327 const SIMINT_DBLTYPE vrr_const_3_over_2pq = SIMINT_MUL(const_3, one_over_2pq);
328 const SIMINT_DBLTYPE vrr_const_4_over_2pq = SIMINT_MUL(const_4, one_over_2pq);
329 const SIMINT_DBLTYPE vrr_const_5_over_2pq = SIMINT_MUL(const_5, one_over_2pq);
330 const SIMINT_DBLTYPE vrr_const_6_over_2pq = SIMINT_MUL(const_6, one_over_2pq);
331 const SIMINT_DBLTYPE vrr_const_7_over_2pq = SIMINT_MUL(const_7, one_over_2pq);
332 const SIMINT_DBLTYPE vrr_const_8_over_2pq = SIMINT_MUL(const_8, one_over_2pq);
333 const SIMINT_DBLTYPE vrr_const_9_over_2pq = SIMINT_MUL(const_9, one_over_2pq);
334 const SIMINT_DBLTYPE vrr_const_10_over_2pq = SIMINT_MUL(const_10, one_over_2pq);
335 const SIMINT_DBLTYPE vrr_const_11_over_2pq = SIMINT_MUL(const_11, one_over_2pq);
336
337
338
339 // Forming PRIM_INT__s_s_p_s[15 * 3];
340 for(n = 0; n < 15; ++n) // loop over orders of auxiliary function
341 {
342
343 PRIM_INT__s_s_p_s[n * 3 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
344 PRIM_INT__s_s_p_s[n * 3 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_p_s[n * 3 + 0]);
345
346 PRIM_INT__s_s_p_s[n * 3 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
347 PRIM_INT__s_s_p_s[n * 3 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_p_s[n * 3 + 1]);
348
349 PRIM_INT__s_s_p_s[n * 3 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
350 PRIM_INT__s_s_p_s[n * 3 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_p_s[n * 3 + 2]);
351
352 }
353
354
355
356 // Forming PRIM_INT__s_s_d_s[14 * 6];
357 for(n = 0; n < 14; ++n) // loop over orders of auxiliary function
358 {
359
360 PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_p_s[n * 3 + 0]);
361 PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_d_s[n * 6 + 0]);
362 PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_s_d_s[n * 6 + 0]);
363
364 PRIM_INT__s_s_d_s[n * 6 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_p_s[n * 3 + 0]);
365 PRIM_INT__s_s_d_s[n * 6 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_d_s[n * 6 + 1]);
366
367 PRIM_INT__s_s_d_s[n * 6 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 0]);
368 PRIM_INT__s_s_d_s[n * 6 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_d_s[n * 6 + 2]);
369
370 PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_p_s[n * 3 + 1]);
371 PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_p_s[(n+1) * 3 + 1], PRIM_INT__s_s_d_s[n * 6 + 3]);
372 PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_s_d_s[n * 6 + 3]);
373
374 PRIM_INT__s_s_d_s[n * 6 + 4] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 1]);
375 PRIM_INT__s_s_d_s[n * 6 + 4] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_p_s[(n+1) * 3 + 1], PRIM_INT__s_s_d_s[n * 6 + 4]);
376
377 PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 2]);
378 PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_p_s[(n+1) * 3 + 2], PRIM_INT__s_s_d_s[n * 6 + 5]);
379 PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_s_d_s[n * 6 + 5]);
380
381 }
382
383
384
385 // Forming PRIM_INT__s_s_f_s[13 * 10];
386 for(n = 0; n < 13; ++n) // loop over orders of auxiliary function
387 {
388
389 PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 0]);
390 PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__s_s_f_s[n * 10 + 0]);
391 PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_FMADD( vrr_const_2_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_p_s[n * 3 + 0]), PRIM_INT__s_s_f_s[n * 10 + 0]);
392
393 PRIM_INT__s_s_f_s[n * 10 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 0]);
394 PRIM_INT__s_s_f_s[n * 10 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__s_s_f_s[n * 10 + 1]);
395
396 PRIM_INT__s_s_f_s[n * 10 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 0]);
397 PRIM_INT__s_s_f_s[n * 10 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__s_s_f_s[n * 10 + 2]);
398
399 PRIM_INT__s_s_f_s[n * 10 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 3]);
400 PRIM_INT__s_s_f_s[n * 10 + 3] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__s_s_f_s[n * 10 + 3]);
401
402 PRIM_INT__s_s_f_s[n * 10 + 4] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 1]);
403 PRIM_INT__s_s_f_s[n * 10 + 4] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 1], PRIM_INT__s_s_f_s[n * 10 + 4]);
404
405 PRIM_INT__s_s_f_s[n * 10 + 5] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 5]);
406 PRIM_INT__s_s_f_s[n * 10 + 5] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__s_s_f_s[n * 10 + 5]);
407
408 PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 3]);
409 PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__s_s_f_s[n * 10 + 6]);
410 PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_FMADD( vrr_const_2_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_p_s[(n+1) * 3 + 1], PRIM_INT__s_s_p_s[n * 3 + 1]), PRIM_INT__s_s_f_s[n * 10 + 6]);
411
412 PRIM_INT__s_s_f_s[n * 10 + 7] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 3]);
413 PRIM_INT__s_s_f_s[n * 10 + 7] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__s_s_f_s[n * 10 + 7]);
414
415 PRIM_INT__s_s_f_s[n * 10 + 8] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 5]);
416 PRIM_INT__s_s_f_s[n * 10 + 8] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__s_s_f_s[n * 10 + 8]);
417
418 PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 5]);
419 PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__s_s_f_s[n * 10 + 9]);
420 PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_FMADD( vrr_const_2_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_p_s[(n+1) * 3 + 2], PRIM_INT__s_s_p_s[n * 3 + 2]), PRIM_INT__s_s_f_s[n * 10 + 9]);
421
422 }
423
424
425 VRR_K_s_s_g_s(
426 PRIM_INT__s_s_g_s,
427 PRIM_INT__s_s_f_s,
428 PRIM_INT__s_s_d_s,
429 Q_PA,
430 a_over_q,
431 aoq_PQ,
432 one_over_2q,
433 12);
434
435
436 VRR_K_s_s_h_s(
437 PRIM_INT__s_s_h_s,
438 PRIM_INT__s_s_g_s,
439 PRIM_INT__s_s_f_s,
440 Q_PA,
441 a_over_q,
442 aoq_PQ,
443 one_over_2q,
444 11);
445
446
447 ostei_general_vrr1_K(6, 10,
448 one_over_2q, a_over_q, aoq_PQ, Q_PA,
449 PRIM_INT__s_s_h_s, PRIM_INT__s_s_g_s, PRIM_INT__s_s_i_s);
450
451
452 ostei_general_vrr_I(1, 0, 6, 0, 4,
453 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
454 PRIM_INT__s_s_i_s, NULL, NULL, PRIM_INT__s_s_h_s, NULL, PRIM_INT__p_s_i_s);
455
456
457 ostei_general_vrr_I(1, 0, 5, 0, 4,
458 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
459 PRIM_INT__s_s_h_s, NULL, NULL, PRIM_INT__s_s_g_s, NULL, PRIM_INT__p_s_h_s);
460
461
462 ostei_general_vrr_I(2, 0, 6, 0, 3,
463 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
464 PRIM_INT__p_s_i_s, PRIM_INT__s_s_i_s, NULL, PRIM_INT__p_s_h_s, NULL, PRIM_INT__d_s_i_s);
465
466
467 VRR_I_p_s_g_s(
468 PRIM_INT__p_s_g_s,
469 PRIM_INT__s_s_g_s,
470 PRIM_INT__s_s_f_s,
471 P_PA,
472 aop_PQ,
473 one_over_2pq,
474 4);
475
476
477 ostei_general_vrr_I(2, 0, 5, 0, 3,
478 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
479 PRIM_INT__p_s_h_s, PRIM_INT__s_s_h_s, NULL, PRIM_INT__p_s_g_s, NULL, PRIM_INT__d_s_h_s);
480
481
482 ostei_general_vrr_I(3, 0, 6, 0, 2,
483 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
484 PRIM_INT__d_s_i_s, PRIM_INT__p_s_i_s, NULL, PRIM_INT__d_s_h_s, NULL, PRIM_INT__f_s_i_s);
485
486
487 ostei_general_vrr1_K(7, 9,
488 one_over_2q, a_over_q, aoq_PQ, Q_PA,
489 PRIM_INT__s_s_i_s, PRIM_INT__s_s_h_s, PRIM_INT__s_s_k_s);
490
491
492 ostei_general_vrr_I(1, 0, 7, 0, 4,
493 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
494 PRIM_INT__s_s_k_s, NULL, NULL, PRIM_INT__s_s_i_s, NULL, PRIM_INT__p_s_k_s);
495
496
497 ostei_general_vrr_I(2, 0, 7, 0, 3,
498 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
499 PRIM_INT__p_s_k_s, PRIM_INT__s_s_k_s, NULL, PRIM_INT__p_s_i_s, NULL, PRIM_INT__d_s_k_s);
500
501
502 ostei_general_vrr_I(3, 0, 7, 0, 2,
503 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
504 PRIM_INT__d_s_k_s, PRIM_INT__p_s_k_s, NULL, PRIM_INT__d_s_i_s, NULL, PRIM_INT__f_s_k_s);
505
506
507 VRR_I_p_s_f_s(
508 PRIM_INT__p_s_f_s,
509 PRIM_INT__s_s_f_s,
510 PRIM_INT__s_s_d_s,
511 P_PA,
512 aop_PQ,
513 one_over_2pq,
514 4);
515
516
517 ostei_general_vrr_I(2, 0, 4, 0, 3,
518 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
519 PRIM_INT__p_s_g_s, PRIM_INT__s_s_g_s, NULL, PRIM_INT__p_s_f_s, NULL, PRIM_INT__d_s_g_s);
520
521
522 ostei_general_vrr_I(3, 0, 5, 0, 2,
523 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
524 PRIM_INT__d_s_h_s, PRIM_INT__p_s_h_s, NULL, PRIM_INT__d_s_g_s, NULL, PRIM_INT__f_s_h_s);
525
526
527 ostei_general_vrr_I(4, 0, 6, 0, 1,
528 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
529 PRIM_INT__f_s_i_s, PRIM_INT__d_s_i_s, NULL, PRIM_INT__f_s_h_s, NULL, PRIM_INT__g_s_i_s);
530
531
532 ostei_general_vrr1_K(8, 8,
533 one_over_2q, a_over_q, aoq_PQ, Q_PA,
534 PRIM_INT__s_s_k_s, PRIM_INT__s_s_i_s, PRIM_INT__s_s_l_s);
535
536
537 ostei_general_vrr_I(1, 0, 8, 0, 4,
538 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
539 PRIM_INT__s_s_l_s, NULL, NULL, PRIM_INT__s_s_k_s, NULL, PRIM_INT__p_s_l_s);
540
541
542 ostei_general_vrr_I(2, 0, 8, 0, 3,
543 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
544 PRIM_INT__p_s_l_s, PRIM_INT__s_s_l_s, NULL, PRIM_INT__p_s_k_s, NULL, PRIM_INT__d_s_l_s);
545
546
547 ostei_general_vrr_I(3, 0, 8, 0, 2,
548 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
549 PRIM_INT__d_s_l_s, PRIM_INT__p_s_l_s, NULL, PRIM_INT__d_s_k_s, NULL, PRIM_INT__f_s_l_s);
550
551
552 ostei_general_vrr_I(4, 0, 7, 0, 1,
553 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
554 PRIM_INT__f_s_k_s, PRIM_INT__d_s_k_s, NULL, PRIM_INT__f_s_i_s, NULL, PRIM_INT__g_s_k_s);
555
556
557 ostei_general_vrr1_K(9, 7,
558 one_over_2q, a_over_q, aoq_PQ, Q_PA,
559 PRIM_INT__s_s_l_s, PRIM_INT__s_s_k_s, PRIM_INT__s_s_m_s);
560
561
562 ostei_general_vrr_I(1, 0, 9, 0, 4,
563 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
564 PRIM_INT__s_s_m_s, NULL, NULL, PRIM_INT__s_s_l_s, NULL, PRIM_INT__p_s_m_s);
565
566
567 ostei_general_vrr_I(2, 0, 9, 0, 3,
568 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
569 PRIM_INT__p_s_m_s, PRIM_INT__s_s_m_s, NULL, PRIM_INT__p_s_l_s, NULL, PRIM_INT__d_s_m_s);
570
571
572 ostei_general_vrr_I(3, 0, 9, 0, 2,
573 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
574 PRIM_INT__d_s_m_s, PRIM_INT__p_s_m_s, NULL, PRIM_INT__d_s_l_s, NULL, PRIM_INT__f_s_m_s);
575
576
577 ostei_general_vrr_I(4, 0, 8, 0, 1,
578 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
579 PRIM_INT__f_s_l_s, PRIM_INT__d_s_l_s, NULL, PRIM_INT__f_s_k_s, NULL, PRIM_INT__g_s_l_s);
580
581
582 ostei_general_vrr1_K(10, 6,
583 one_over_2q, a_over_q, aoq_PQ, Q_PA,
584 PRIM_INT__s_s_m_s, PRIM_INT__s_s_l_s, PRIM_INT__s_s_n_s);
585
586
587 ostei_general_vrr_I(1, 0, 10, 0, 4,
588 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
589 PRIM_INT__s_s_n_s, NULL, NULL, PRIM_INT__s_s_m_s, NULL, PRIM_INT__p_s_n_s);
590
591
592 ostei_general_vrr_I(2, 0, 10, 0, 3,
593 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
594 PRIM_INT__p_s_n_s, PRIM_INT__s_s_n_s, NULL, PRIM_INT__p_s_m_s, NULL, PRIM_INT__d_s_n_s);
595
596
597 ostei_general_vrr_I(3, 0, 10, 0, 2,
598 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
599 PRIM_INT__d_s_n_s, PRIM_INT__p_s_n_s, NULL, PRIM_INT__d_s_m_s, NULL, PRIM_INT__f_s_n_s);
600
601
602 ostei_general_vrr_I(4, 0, 9, 0, 1,
603 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
604 PRIM_INT__f_s_m_s, PRIM_INT__d_s_m_s, NULL, PRIM_INT__f_s_l_s, NULL, PRIM_INT__g_s_m_s);
605
606
607 ostei_general_vrr1_K(11, 5,
608 one_over_2q, a_over_q, aoq_PQ, Q_PA,
609 PRIM_INT__s_s_n_s, PRIM_INT__s_s_m_s, PRIM_INT__s_s_o_s);
610
611
612 ostei_general_vrr_I(1, 0, 11, 0, 4,
613 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
614 PRIM_INT__s_s_o_s, NULL, NULL, PRIM_INT__s_s_n_s, NULL, PRIM_INT__p_s_o_s);
615
616
617 ostei_general_vrr_I(2, 0, 11, 0, 3,
618 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
619 PRIM_INT__p_s_o_s, PRIM_INT__s_s_o_s, NULL, PRIM_INT__p_s_n_s, NULL, PRIM_INT__d_s_o_s);
620
621
622 ostei_general_vrr_I(3, 0, 11, 0, 2,
623 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
624 PRIM_INT__d_s_o_s, PRIM_INT__p_s_o_s, NULL, PRIM_INT__d_s_n_s, NULL, PRIM_INT__f_s_o_s);
625
626
627 ostei_general_vrr_I(4, 0, 10, 0, 1,
628 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
629 PRIM_INT__f_s_n_s, PRIM_INT__d_s_n_s, NULL, PRIM_INT__f_s_m_s, NULL, PRIM_INT__g_s_n_s);
630
631
632 ostei_general_vrr_I(4, 0, 11, 0, 1,
633 one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
634 PRIM_INT__f_s_o_s, PRIM_INT__d_s_o_s, NULL, PRIM_INT__f_s_n_s, NULL, PRIM_INT__g_s_o_s);
635
636
637
638
639 ////////////////////////////////////
640 // Accumulate contracted integrals
641 ////////////////////////////////////
642 if(lastoffset == 0)
643 {
644 contract_all(280, PRIM_INT__f_s_i_s, PRIM_PTR_INT__f_s_i_s);
645 contract_all(360, PRIM_INT__f_s_k_s, PRIM_PTR_INT__f_s_k_s);
646 contract_all(450, PRIM_INT__f_s_l_s, PRIM_PTR_INT__f_s_l_s);
647 contract_all(550, PRIM_INT__f_s_m_s, PRIM_PTR_INT__f_s_m_s);
648 contract_all(660, PRIM_INT__f_s_n_s, PRIM_PTR_INT__f_s_n_s);
649 contract_all(780, PRIM_INT__f_s_o_s, PRIM_PTR_INT__f_s_o_s);
650 contract_all(420, PRIM_INT__g_s_i_s, PRIM_PTR_INT__g_s_i_s);
651 contract_all(540, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
652 contract_all(675, PRIM_INT__g_s_l_s, PRIM_PTR_INT__g_s_l_s);
653 contract_all(825, PRIM_INT__g_s_m_s, PRIM_PTR_INT__g_s_m_s);
654 contract_all(990, PRIM_INT__g_s_n_s, PRIM_PTR_INT__g_s_n_s);
655 contract_all(1170, PRIM_INT__g_s_o_s, PRIM_PTR_INT__g_s_o_s);
656 }
657 else
658 {
659 contract(280, shelloffsets, PRIM_INT__f_s_i_s, PRIM_PTR_INT__f_s_i_s);
660 contract(360, shelloffsets, PRIM_INT__f_s_k_s, PRIM_PTR_INT__f_s_k_s);
661 contract(450, shelloffsets, PRIM_INT__f_s_l_s, PRIM_PTR_INT__f_s_l_s);
662 contract(550, shelloffsets, PRIM_INT__f_s_m_s, PRIM_PTR_INT__f_s_m_s);
663 contract(660, shelloffsets, PRIM_INT__f_s_n_s, PRIM_PTR_INT__f_s_n_s);
664 contract(780, shelloffsets, PRIM_INT__f_s_o_s, PRIM_PTR_INT__f_s_o_s);
665 contract(420, shelloffsets, PRIM_INT__g_s_i_s, PRIM_PTR_INT__g_s_i_s);
666 contract(540, shelloffsets, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
667 contract(675, shelloffsets, PRIM_INT__g_s_l_s, PRIM_PTR_INT__g_s_l_s);
668 contract(825, shelloffsets, PRIM_INT__g_s_m_s, PRIM_PTR_INT__g_s_m_s);
669 contract(990, shelloffsets, PRIM_INT__g_s_n_s, PRIM_PTR_INT__g_s_n_s);
670 contract(1170, shelloffsets, PRIM_INT__g_s_o_s, PRIM_PTR_INT__g_s_o_s);
671 PRIM_PTR_INT__f_s_i_s += lastoffset*280;
672 PRIM_PTR_INT__f_s_k_s += lastoffset*360;
673 PRIM_PTR_INT__f_s_l_s += lastoffset*450;
674 PRIM_PTR_INT__f_s_m_s += lastoffset*550;
675 PRIM_PTR_INT__f_s_n_s += lastoffset*660;
676 PRIM_PTR_INT__f_s_o_s += lastoffset*780;
677 PRIM_PTR_INT__g_s_i_s += lastoffset*420;
678 PRIM_PTR_INT__g_s_k_s += lastoffset*540;
679 PRIM_PTR_INT__g_s_l_s += lastoffset*675;
680 PRIM_PTR_INT__g_s_m_s += lastoffset*825;
681 PRIM_PTR_INT__g_s_n_s += lastoffset*990;
682 PRIM_PTR_INT__g_s_o_s += lastoffset*1170;
683 }
684
685 } // close loop over j
686 } // close loop over i
687
688 //Advance to the next batch
689 jstart = SIMINT_SIMD_ROUND(jend);
690
691 //////////////////////////////////////////////
692 // Contracted integrals: Horizontal recurrance
693 //////////////////////////////////////////////
694
695
696 const double hAB[3] = { P.AB_x[ab], P.AB_y[ab], P.AB_z[ab] };
697
698
699 for(abcd = 0; abcd < nshellbatch; ++abcd, ++real_abcd)
700 {
701 const double hCD[3] = { Q.AB_x[cd+abcd], Q.AB_y[cd+abcd], Q.AB_z[cd+abcd] };
702
703 // set up HRR pointers
704 double const * restrict HRR_INT__f_s_i_s = INT__f_s_i_s + abcd * 280;
705 double const * restrict HRR_INT__f_s_k_s = INT__f_s_k_s + abcd * 360;
706 double const * restrict HRR_INT__f_s_l_s = INT__f_s_l_s + abcd * 450;
707 double const * restrict HRR_INT__f_s_m_s = INT__f_s_m_s + abcd * 550;
708 double const * restrict HRR_INT__f_s_n_s = INT__f_s_n_s + abcd * 660;
709 double const * restrict HRR_INT__f_s_o_s = INT__f_s_o_s + abcd * 780;
710 double const * restrict HRR_INT__g_s_i_s = INT__g_s_i_s + abcd * 420;
711 double const * restrict HRR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
712 double const * restrict HRR_INT__g_s_l_s = INT__g_s_l_s + abcd * 675;
713 double const * restrict HRR_INT__g_s_m_s = INT__g_s_m_s + abcd * 825;
714 double const * restrict HRR_INT__g_s_n_s = INT__g_s_n_s + abcd * 990;
715 double const * restrict HRR_INT__g_s_o_s = INT__g_s_o_s + abcd * 1170;
716 double * restrict HRR_INT__f_p_i_h = INT__f_p_i_h + real_abcd * 17640;
717
718 // form INT__f_p_i_s
719 HRR_J_f_p(
720 HRR_INT__f_p_i_s,
721 HRR_INT__f_s_i_s,
722 HRR_INT__g_s_i_s,
723 hAB, 28);
724
725 // form INT__f_p_k_s
726 HRR_J_f_p(
727 HRR_INT__f_p_k_s,
728 HRR_INT__f_s_k_s,
729 HRR_INT__g_s_k_s,
730 hAB, 36);
731
732 // form INT__f_p_l_s
733 HRR_J_f_p(
734 HRR_INT__f_p_l_s,
735 HRR_INT__f_s_l_s,
736 HRR_INT__g_s_l_s,
737 hAB, 45);
738
739 // form INT__f_p_m_s
740 HRR_J_f_p(
741 HRR_INT__f_p_m_s,
742 HRR_INT__f_s_m_s,
743 HRR_INT__g_s_m_s,
744 hAB, 55);
745
746 // form INT__f_p_n_s
747 HRR_J_f_p(
748 HRR_INT__f_p_n_s,
749 HRR_INT__f_s_n_s,
750 HRR_INT__g_s_n_s,
751 hAB, 66);
752
753 // form INT__f_p_o_s
754 HRR_J_f_p(
755 HRR_INT__f_p_o_s,
756 HRR_INT__f_s_o_s,
757 HRR_INT__g_s_o_s,
758 hAB, 78);
759
760 // form INT__f_p_i_p
761 ostei_general_hrr_L(3, 1, 6, 1, hCD, HRR_INT__f_p_k_s, HRR_INT__f_p_i_s, HRR_INT__f_p_i_p);
762
763 // form INT__f_p_k_p
764 ostei_general_hrr_L(3, 1, 7, 1, hCD, HRR_INT__f_p_l_s, HRR_INT__f_p_k_s, HRR_INT__f_p_k_p);
765
766 // form INT__f_p_l_p
767 ostei_general_hrr_L(3, 1, 8, 1, hCD, HRR_INT__f_p_m_s, HRR_INT__f_p_l_s, HRR_INT__f_p_l_p);
768
769 // form INT__f_p_m_p
770 ostei_general_hrr_L(3, 1, 9, 1, hCD, HRR_INT__f_p_n_s, HRR_INT__f_p_m_s, HRR_INT__f_p_m_p);
771
772 // form INT__f_p_n_p
773 ostei_general_hrr_L(3, 1, 10, 1, hCD, HRR_INT__f_p_o_s, HRR_INT__f_p_n_s, HRR_INT__f_p_n_p);
774
775 // form INT__f_p_i_d
776 ostei_general_hrr_L(3, 1, 6, 2, hCD, HRR_INT__f_p_k_p, HRR_INT__f_p_i_p, HRR_INT__f_p_i_d);
777
778 // form INT__f_p_k_d
779 ostei_general_hrr_L(3, 1, 7, 2, hCD, HRR_INT__f_p_l_p, HRR_INT__f_p_k_p, HRR_INT__f_p_k_d);
780
781 // form INT__f_p_l_d
782 ostei_general_hrr_L(3, 1, 8, 2, hCD, HRR_INT__f_p_m_p, HRR_INT__f_p_l_p, HRR_INT__f_p_l_d);
783
784 // form INT__f_p_m_d
785 ostei_general_hrr_L(3, 1, 9, 2, hCD, HRR_INT__f_p_n_p, HRR_INT__f_p_m_p, HRR_INT__f_p_m_d);
786
787 // form INT__f_p_i_f
788 ostei_general_hrr_L(3, 1, 6, 3, hCD, HRR_INT__f_p_k_d, HRR_INT__f_p_i_d, HRR_INT__f_p_i_f);
789
790 // form INT__f_p_k_f
791 ostei_general_hrr_L(3, 1, 7, 3, hCD, HRR_INT__f_p_l_d, HRR_INT__f_p_k_d, HRR_INT__f_p_k_f);
792
793 // form INT__f_p_l_f
794 ostei_general_hrr_L(3, 1, 8, 3, hCD, HRR_INT__f_p_m_d, HRR_INT__f_p_l_d, HRR_INT__f_p_l_f);
795
796 // form INT__f_p_i_g
797 ostei_general_hrr_L(3, 1, 6, 4, hCD, HRR_INT__f_p_k_f, HRR_INT__f_p_i_f, HRR_INT__f_p_i_g);
798
799 // form INT__f_p_k_g
800 ostei_general_hrr_L(3, 1, 7, 4, hCD, HRR_INT__f_p_l_f, HRR_INT__f_p_k_f, HRR_INT__f_p_k_g);
801
802 // form INT__f_p_i_h
803 ostei_general_hrr_L(3, 1, 6, 5, hCD, HRR_INT__f_p_k_g, HRR_INT__f_p_i_g, HRR_INT__f_p_i_h);
804
805
806 } // close HRR loop
807
808
809 } // close loop cdbatch
810
811 istart = iend;
812 } // close loop over ab
813
814 return P.nshell12_clip * Q.nshell12_clip;
815 }
816
ostei_p_f_i_h(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__p_f_i_h)817 int ostei_p_f_i_h(struct simint_multi_shellpair const P,
818 struct simint_multi_shellpair const Q,
819 double screen_tol,
820 double * const restrict work,
821 double * const restrict INT__p_f_i_h)
822 {
823 double P_AB[3*P.nshell12];
824 struct simint_multi_shellpair P_tmp = P;
825 P_tmp.PA_x = P.PB_x; P_tmp.PA_y = P.PB_y; P_tmp.PA_z = P.PB_z;
826 P_tmp.PB_x = P.PA_x; P_tmp.PB_y = P.PA_y; P_tmp.PB_z = P.PA_z;
827 P_tmp.AB_x = P_AB;
828 P_tmp.AB_y = P_AB + P.nshell12;
829 P_tmp.AB_z = P_AB + 2*P.nshell12;
830
831 for(int i = 0; i < P.nshell12; i++)
832 {
833 P_tmp.AB_x[i] = -P.AB_x[i];
834 P_tmp.AB_y[i] = -P.AB_y[i];
835 P_tmp.AB_z[i] = -P.AB_z[i];
836 }
837
838 int ret = ostei_f_p_i_h(P_tmp, Q, screen_tol, work, INT__p_f_i_h);
839 double buffer[17640] SIMINT_ALIGN_ARRAY_DBL;
840
841 for(int q = 0; q < ret; q++)
842 {
843 int idx = 0;
844 for(int a = 0; a < 3; ++a)
845 for(int b = 0; b < 10; ++b)
846 for(int c = 0; c < 28; ++c)
847 for(int d = 0; d < 21; ++d)
848 buffer[idx++] = INT__p_f_i_h[q*17640+b*1764+a*588+c*21+d];
849
850 memcpy(INT__p_f_i_h+q*17640, buffer, 17640*sizeof(double));
851 }
852
853 return ret;
854 }
855
ostei_f_p_h_i(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__f_p_h_i)856 int ostei_f_p_h_i(struct simint_multi_shellpair const P,
857 struct simint_multi_shellpair const Q,
858 double screen_tol,
859 double * const restrict work,
860 double * const restrict INT__f_p_h_i)
861 {
862 double Q_AB[3*Q.nshell12];
863 struct simint_multi_shellpair Q_tmp = Q;
864 Q_tmp.PA_x = Q.PB_x; Q_tmp.PA_y = Q.PB_y; Q_tmp.PA_z = Q.PB_z;
865 Q_tmp.PB_x = Q.PA_x; Q_tmp.PB_y = Q.PA_y; Q_tmp.PB_z = Q.PA_z;
866 Q_tmp.AB_x = Q_AB;
867 Q_tmp.AB_y = Q_AB + Q.nshell12;
868 Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
869
870 for(int i = 0; i < Q.nshell12; i++)
871 {
872 Q_tmp.AB_x[i] = -Q.AB_x[i];
873 Q_tmp.AB_y[i] = -Q.AB_y[i];
874 Q_tmp.AB_z[i] = -Q.AB_z[i];
875 }
876
877 int ret = ostei_f_p_i_h(P, Q_tmp, screen_tol, work, INT__f_p_h_i);
878 double buffer[17640] SIMINT_ALIGN_ARRAY_DBL;
879
880 for(int q = 0; q < ret; q++)
881 {
882 int idx = 0;
883 for(int a = 0; a < 10; ++a)
884 for(int b = 0; b < 3; ++b)
885 for(int c = 0; c < 21; ++c)
886 for(int d = 0; d < 28; ++d)
887 buffer[idx++] = INT__f_p_h_i[q*17640+a*1764+b*588+d*21+c];
888
889 memcpy(INT__f_p_h_i+q*17640, buffer, 17640*sizeof(double));
890 }
891
892 return ret;
893 }
894
ostei_p_f_h_i(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__p_f_h_i)895 int ostei_p_f_h_i(struct simint_multi_shellpair const P,
896 struct simint_multi_shellpair const Q,
897 double screen_tol,
898 double * const restrict work,
899 double * const restrict INT__p_f_h_i)
900 {
901 double P_AB[3*P.nshell12];
902 struct simint_multi_shellpair P_tmp = P;
903 P_tmp.PA_x = P.PB_x; P_tmp.PA_y = P.PB_y; P_tmp.PA_z = P.PB_z;
904 P_tmp.PB_x = P.PA_x; P_tmp.PB_y = P.PA_y; P_tmp.PB_z = P.PA_z;
905 P_tmp.AB_x = P_AB;
906 P_tmp.AB_y = P_AB + P.nshell12;
907 P_tmp.AB_z = P_AB + 2*P.nshell12;
908
909 for(int i = 0; i < P.nshell12; i++)
910 {
911 P_tmp.AB_x[i] = -P.AB_x[i];
912 P_tmp.AB_y[i] = -P.AB_y[i];
913 P_tmp.AB_z[i] = -P.AB_z[i];
914 }
915
916 double Q_AB[3*Q.nshell12];
917 struct simint_multi_shellpair Q_tmp = Q;
918 Q_tmp.PA_x = Q.PB_x; Q_tmp.PA_y = Q.PB_y; Q_tmp.PA_z = Q.PB_z;
919 Q_tmp.PB_x = Q.PA_x; Q_tmp.PB_y = Q.PA_y; Q_tmp.PB_z = Q.PA_z;
920 Q_tmp.AB_x = Q_AB;
921 Q_tmp.AB_y = Q_AB + Q.nshell12;
922 Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
923
924 for(int i = 0; i < Q.nshell12; i++)
925 {
926 Q_tmp.AB_x[i] = -Q.AB_x[i];
927 Q_tmp.AB_y[i] = -Q.AB_y[i];
928 Q_tmp.AB_z[i] = -Q.AB_z[i];
929 }
930
931 int ret = ostei_f_p_i_h(P_tmp, Q_tmp, screen_tol, work, INT__p_f_h_i);
932 double buffer[17640] SIMINT_ALIGN_ARRAY_DBL;
933
934 for(int q = 0; q < ret; q++)
935 {
936 int idx = 0;
937 for(int a = 0; a < 3; ++a)
938 for(int b = 0; b < 10; ++b)
939 for(int c = 0; c < 21; ++c)
940 for(int d = 0; d < 28; ++d)
941 buffer[idx++] = INT__p_f_h_i[q*17640+b*1764+a*588+d*21+c];
942
943 memcpy(INT__p_f_h_i+q*17640, buffer, 17640*sizeof(double));
944 }
945
946 return ret;
947 }
948
949