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_d_f_p_f(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__d_f_p_f)8 int ostei_d_f_p_f(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__d_f_p_f)
13 {
14
15 SIMINT_ASSUME_ALIGN_DBL(work);
16 SIMINT_ASSUME_ALIGN_DBL(INT__d_f_p_f);
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__s_f_s_f = work + (SIMINT_NSHELL_SIMD * 0);
30 double * const INT__s_f_s_g = work + (SIMINT_NSHELL_SIMD * 100);
31 double * const INT__s_g_s_f = work + (SIMINT_NSHELL_SIMD * 250);
32 double * const INT__s_g_s_g = work + (SIMINT_NSHELL_SIMD * 400);
33 double * const INT__s_h_s_f = work + (SIMINT_NSHELL_SIMD * 625);
34 double * const INT__s_h_s_g = work + (SIMINT_NSHELL_SIMD * 835);
35 SIMINT_DBLTYPE * const primwork = (SIMINT_DBLTYPE *)(work + SIMINT_NSHELL_SIMD*1150);
36 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_s = primwork + 0;
37 SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_p = primwork + 10;
38 SIMINT_DBLTYPE * const restrict PRIM_INT__s_p_s_s = primwork + 22;
39 SIMINT_DBLTYPE * const restrict PRIM_INT__s_p_s_p = primwork + 49;
40 SIMINT_DBLTYPE * const restrict PRIM_INT__s_p_s_d = primwork + 85;
41 SIMINT_DBLTYPE * const restrict PRIM_INT__s_d_s_s = primwork + 139;
42 SIMINT_DBLTYPE * const restrict PRIM_INT__s_d_s_p = primwork + 187;
43 SIMINT_DBLTYPE * const restrict PRIM_INT__s_d_s_d = primwork + 259;
44 SIMINT_DBLTYPE * const restrict PRIM_INT__s_d_s_f = primwork + 367;
45 SIMINT_DBLTYPE * const restrict PRIM_INT__s_f_s_s = primwork + 487;
46 SIMINT_DBLTYPE * const restrict PRIM_INT__s_f_s_p = primwork + 557;
47 SIMINT_DBLTYPE * const restrict PRIM_INT__s_f_s_d = primwork + 677;
48 SIMINT_DBLTYPE * const restrict PRIM_INT__s_f_s_f = primwork + 857;
49 SIMINT_DBLTYPE * const restrict PRIM_INT__s_f_s_g = primwork + 1057;
50 SIMINT_DBLTYPE * const restrict PRIM_INT__s_g_s_s = primwork + 1207;
51 SIMINT_DBLTYPE * const restrict PRIM_INT__s_g_s_p = primwork + 1297;
52 SIMINT_DBLTYPE * const restrict PRIM_INT__s_g_s_d = primwork + 1477;
53 SIMINT_DBLTYPE * const restrict PRIM_INT__s_g_s_f = primwork + 1747;
54 SIMINT_DBLTYPE * const restrict PRIM_INT__s_g_s_g = primwork + 2047;
55 SIMINT_DBLTYPE * const restrict PRIM_INT__s_h_s_s = primwork + 2272;
56 SIMINT_DBLTYPE * const restrict PRIM_INT__s_h_s_p = primwork + 2377;
57 SIMINT_DBLTYPE * const restrict PRIM_INT__s_h_s_d = primwork + 2629;
58 SIMINT_DBLTYPE * const restrict PRIM_INT__s_h_s_f = primwork + 3007;
59 SIMINT_DBLTYPE * const restrict PRIM_INT__s_h_s_g = primwork + 3427;
60 double * const hrrwork = (double *)(primwork + 3742);
61 double * const HRR_INT__p_f_s_f = hrrwork + 0;
62 double * const HRR_INT__p_f_s_g = hrrwork + 300;
63 double * const HRR_INT__p_g_s_f = hrrwork + 750;
64 double * const HRR_INT__p_g_s_g = hrrwork + 1200;
65 double * const HRR_INT__d_f_s_f = hrrwork + 1875;
66 double * const HRR_INT__d_f_s_g = hrrwork + 2475;
67
68
69 // Create constants
70 const SIMINT_DBLTYPE const_1 = SIMINT_DBLSET1(1);
71 const SIMINT_DBLTYPE const_2 = SIMINT_DBLSET1(2);
72 const SIMINT_DBLTYPE const_3 = SIMINT_DBLSET1(3);
73 const SIMINT_DBLTYPE const_4 = SIMINT_DBLSET1(4);
74 const SIMINT_DBLTYPE const_5 = SIMINT_DBLSET1(5);
75 const SIMINT_DBLTYPE one_half = SIMINT_DBLSET1(0.5);
76
77
78 ////////////////////////////////////////
79 // Loop over shells and primitives
80 ////////////////////////////////////////
81
82 real_abcd = 0;
83 istart = 0;
84 for(ab = 0; ab < P.nshell12_clip; ++ab)
85 {
86 const int iend = istart + P.nprim12[ab];
87
88 cd = 0;
89 jstart = 0;
90
91 for(cd = 0; cd < Q.nshell12_clip; cd += SIMINT_NSHELL_SIMD)
92 {
93 const int nshellbatch = ((cd + SIMINT_NSHELL_SIMD) > Q.nshell12_clip) ? Q.nshell12_clip - cd : SIMINT_NSHELL_SIMD;
94 int jend = jstart;
95 for(i = 0; i < nshellbatch; i++)
96 jend += Q.nprim12[cd+i];
97
98 // Clear the beginning of the workspace (where we are accumulating integrals)
99 memset(work, 0, SIMINT_NSHELL_SIMD * 1150 * sizeof(double));
100 abcd = 0;
101
102
103 for(i = istart; i < iend; ++i)
104 {
105 SIMINT_DBLTYPE bra_screen_max; // only used if check_screen
106
107 if(check_screen)
108 {
109 // Skip this whole thing if always insignificant
110 if((P.screen[i] * Q.screen_max) < screen_tol)
111 continue;
112 bra_screen_max = SIMINT_DBLSET1(P.screen[i]);
113 }
114
115 icd = 0;
116 iprimcd = 0;
117 nprim_icd = Q.nprim12[cd];
118 double * restrict PRIM_PTR_INT__s_f_s_f = INT__s_f_s_f + abcd * 100;
119 double * restrict PRIM_PTR_INT__s_f_s_g = INT__s_f_s_g + abcd * 150;
120 double * restrict PRIM_PTR_INT__s_g_s_f = INT__s_g_s_f + abcd * 150;
121 double * restrict PRIM_PTR_INT__s_g_s_g = INT__s_g_s_g + abcd * 225;
122 double * restrict PRIM_PTR_INT__s_h_s_f = INT__s_h_s_f + abcd * 210;
123 double * restrict PRIM_PTR_INT__s_h_s_g = INT__s_h_s_g + abcd * 315;
124
125
126
127 // Load these one per loop over i
128 const SIMINT_DBLTYPE P_alpha = SIMINT_DBLSET1(P.alpha[i]);
129 const SIMINT_DBLTYPE P_prefac = SIMINT_DBLSET1(P.prefac[i]);
130 const SIMINT_DBLTYPE Pxyz[3] = { SIMINT_DBLSET1(P.x[i]), SIMINT_DBLSET1(P.y[i]), SIMINT_DBLSET1(P.z[i]) };
131
132 const SIMINT_DBLTYPE P_PB[3] = { SIMINT_DBLSET1(P.PB_x[i]), SIMINT_DBLSET1(P.PB_y[i]), SIMINT_DBLSET1(P.PB_z[i]) };
133
134 for(j = jstart; j < jend; j += SIMINT_SIMD_LEN)
135 {
136 // calculate the shell offsets
137 // these are the offset from the shell pointed to by cd
138 // for each element
139 int shelloffsets[SIMINT_SIMD_LEN] = {0};
140 int lastoffset = 0;
141 const int nlane = ( ((j + SIMINT_SIMD_LEN) < jend) ? SIMINT_SIMD_LEN : (jend - j));
142
143 if((iprimcd + SIMINT_SIMD_LEN) >= nprim_icd)
144 {
145 // Handle if the first element of the vector is a new shell
146 if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
147 {
148 nprim_icd += Q.nprim12[cd + (++icd)];
149 PRIM_PTR_INT__s_f_s_f += 100;
150 PRIM_PTR_INT__s_f_s_g += 150;
151 PRIM_PTR_INT__s_g_s_f += 150;
152 PRIM_PTR_INT__s_g_s_g += 225;
153 PRIM_PTR_INT__s_h_s_f += 210;
154 PRIM_PTR_INT__s_h_s_g += 315;
155 }
156 iprimcd++;
157 for(n = 1; n < SIMINT_SIMD_LEN; ++n)
158 {
159 if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
160 {
161 shelloffsets[n] = shelloffsets[n-1] + 1;
162 lastoffset++;
163 nprim_icd += Q.nprim12[cd + (++icd)];
164 }
165 else
166 shelloffsets[n] = shelloffsets[n-1];
167 iprimcd++;
168 }
169 }
170 else
171 iprimcd += SIMINT_SIMD_LEN;
172
173 // Do we have to compute this vector (or has it been screened out)?
174 // (not_screened != 0 means we have to do this vector)
175 if(check_screen)
176 {
177 const double vmax = vector_max(SIMINT_MUL(bra_screen_max, SIMINT_DBLLOAD(Q.screen, j)));
178 if(vmax < screen_tol)
179 {
180 PRIM_PTR_INT__s_f_s_f += lastoffset*100;
181 PRIM_PTR_INT__s_f_s_g += lastoffset*150;
182 PRIM_PTR_INT__s_g_s_f += lastoffset*150;
183 PRIM_PTR_INT__s_g_s_g += lastoffset*225;
184 PRIM_PTR_INT__s_h_s_f += lastoffset*210;
185 PRIM_PTR_INT__s_h_s_g += lastoffset*315;
186 continue;
187 }
188 }
189
190 const SIMINT_DBLTYPE Q_alpha = SIMINT_DBLLOAD(Q.alpha, j);
191 const SIMINT_DBLTYPE PQalpha_mul = SIMINT_MUL(P_alpha, Q_alpha);
192 const SIMINT_DBLTYPE PQalpha_sum = SIMINT_ADD(P_alpha, Q_alpha);
193 const SIMINT_DBLTYPE one_over_PQalpha_sum = SIMINT_DIV(const_1, PQalpha_sum);
194
195
196 /* construct R2 = (Px - Qx)**2 + (Py - Qy)**2 + (Pz -Qz)**2 */
197 SIMINT_DBLTYPE PQ[3];
198 PQ[0] = SIMINT_SUB(Pxyz[0], SIMINT_DBLLOAD(Q.x, j));
199 PQ[1] = SIMINT_SUB(Pxyz[1], SIMINT_DBLLOAD(Q.y, j));
200 PQ[2] = SIMINT_SUB(Pxyz[2], SIMINT_DBLLOAD(Q.z, j));
201 SIMINT_DBLTYPE R2 = SIMINT_MUL(PQ[0], PQ[0]);
202 R2 = SIMINT_FMADD(PQ[1], PQ[1], R2);
203 R2 = SIMINT_FMADD(PQ[2], PQ[2], R2);
204
205 const SIMINT_DBLTYPE alpha = SIMINT_MUL(PQalpha_mul, one_over_PQalpha_sum); // alpha from MEST
206 const SIMINT_DBLTYPE one_over_p = SIMINT_DIV(const_1, P_alpha);
207 const SIMINT_DBLTYPE one_over_q = SIMINT_DIV(const_1, Q_alpha);
208 const SIMINT_DBLTYPE one_over_2p = SIMINT_MUL(one_half, one_over_p);
209 const SIMINT_DBLTYPE one_over_2q = SIMINT_MUL(one_half, one_over_q);
210 const SIMINT_DBLTYPE one_over_2pq = SIMINT_MUL(one_half, one_over_PQalpha_sum);
211 const SIMINT_DBLTYPE Q_PB[3] = { SIMINT_DBLLOAD(Q.PB_x, j), SIMINT_DBLLOAD(Q.PB_y, j), SIMINT_DBLLOAD(Q.PB_z, j) };
212
213 // NOTE: Minus sign!
214 const SIMINT_DBLTYPE a_over_p = SIMINT_MUL(SIMINT_NEG(alpha), one_over_p);
215 SIMINT_DBLTYPE aop_PQ[3];
216 aop_PQ[0] = SIMINT_MUL(a_over_p, PQ[0]);
217 aop_PQ[1] = SIMINT_MUL(a_over_p, PQ[1]);
218 aop_PQ[2] = SIMINT_MUL(a_over_p, PQ[2]);
219
220 SIMINT_DBLTYPE a_over_q = SIMINT_MUL(alpha, one_over_q);
221 SIMINT_DBLTYPE aoq_PQ[3];
222 aoq_PQ[0] = SIMINT_MUL(a_over_q, PQ[0]);
223 aoq_PQ[1] = SIMINT_MUL(a_over_q, PQ[1]);
224 aoq_PQ[2] = SIMINT_MUL(a_over_q, PQ[2]);
225 // Put a minus sign here so we don't have to in RR routines
226 a_over_q = SIMINT_NEG(a_over_q);
227
228
229 //////////////////////////////////////////////
230 // Fjt function section
231 // Maximum v value: 9
232 //////////////////////////////////////////////
233 // The parameter to the Fjt function
234 const SIMINT_DBLTYPE F_x = SIMINT_MUL(R2, alpha);
235
236
237 const SIMINT_DBLTYPE Q_prefac = mask_load(nlane, Q.prefac + j);
238
239
240 boys_F_split(PRIM_INT__s_s_s_s, F_x, 9);
241 SIMINT_DBLTYPE prefac = SIMINT_SQRT(one_over_PQalpha_sum);
242 prefac = SIMINT_MUL(SIMINT_MUL(P_prefac, Q_prefac), prefac);
243 for(n = 0; n <= 9; n++)
244 PRIM_INT__s_s_s_s[n] = SIMINT_MUL(PRIM_INT__s_s_s_s[n], prefac);
245
246 //////////////////////////////////////////////
247 // Primitive integrals: Vertical recurrance
248 //////////////////////////////////////////////
249
250 const SIMINT_DBLTYPE vrr_const_1_over_2p = one_over_2p;
251 const SIMINT_DBLTYPE vrr_const_2_over_2p = SIMINT_MUL(const_2, one_over_2p);
252 const SIMINT_DBLTYPE vrr_const_3_over_2p = SIMINT_MUL(const_3, one_over_2p);
253 const SIMINT_DBLTYPE vrr_const_4_over_2p = SIMINT_MUL(const_4, one_over_2p);
254 const SIMINT_DBLTYPE vrr_const_1_over_2q = one_over_2q;
255 const SIMINT_DBLTYPE vrr_const_2_over_2q = SIMINT_MUL(const_2, one_over_2q);
256 const SIMINT_DBLTYPE vrr_const_3_over_2q = SIMINT_MUL(const_3, one_over_2q);
257 const SIMINT_DBLTYPE vrr_const_1_over_2pq = one_over_2pq;
258 const SIMINT_DBLTYPE vrr_const_2_over_2pq = SIMINT_MUL(const_2, one_over_2pq);
259 const SIMINT_DBLTYPE vrr_const_3_over_2pq = SIMINT_MUL(const_3, one_over_2pq);
260 const SIMINT_DBLTYPE vrr_const_4_over_2pq = SIMINT_MUL(const_4, one_over_2pq);
261 const SIMINT_DBLTYPE vrr_const_5_over_2pq = SIMINT_MUL(const_5, one_over_2pq);
262
263
264
265 // Forming PRIM_INT__s_p_s_s[9 * 3];
266 for(n = 0; n < 9; ++n) // loop over orders of auxiliary function
267 {
268
269 PRIM_INT__s_p_s_s[n * 3 + 0] = SIMINT_MUL(P_PB[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
270 PRIM_INT__s_p_s_s[n * 3 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_p_s_s[n * 3 + 0]);
271
272 PRIM_INT__s_p_s_s[n * 3 + 1] = SIMINT_MUL(P_PB[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
273 PRIM_INT__s_p_s_s[n * 3 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_p_s_s[n * 3 + 1]);
274
275 PRIM_INT__s_p_s_s[n * 3 + 2] = SIMINT_MUL(P_PB[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
276 PRIM_INT__s_p_s_s[n * 3 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_p_s_s[n * 3 + 2]);
277
278 }
279
280
281
282 // Forming PRIM_INT__s_d_s_s[8 * 6];
283 for(n = 0; n < 8; ++n) // loop over orders of auxiliary function
284 {
285
286 PRIM_INT__s_d_s_s[n * 6 + 0] = SIMINT_MUL(P_PB[0], PRIM_INT__s_p_s_s[n * 3 + 0]);
287 PRIM_INT__s_d_s_s[n * 6 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_d_s_s[n * 6 + 0]);
288 PRIM_INT__s_d_s_s[n * 6 + 0] = SIMINT_FMADD( vrr_const_1_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_d_s_s[n * 6 + 0]);
289
290 PRIM_INT__s_d_s_s[n * 6 + 1] = SIMINT_MUL(P_PB[1], PRIM_INT__s_p_s_s[n * 3 + 0]);
291 PRIM_INT__s_d_s_s[n * 6 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_d_s_s[n * 6 + 1]);
292
293 PRIM_INT__s_d_s_s[n * 6 + 2] = SIMINT_MUL(P_PB[2], PRIM_INT__s_p_s_s[n * 3 + 0]);
294 PRIM_INT__s_d_s_s[n * 6 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_d_s_s[n * 6 + 2]);
295
296 PRIM_INT__s_d_s_s[n * 6 + 3] = SIMINT_MUL(P_PB[1], PRIM_INT__s_p_s_s[n * 3 + 1]);
297 PRIM_INT__s_d_s_s[n * 6 + 3] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_d_s_s[n * 6 + 3]);
298 PRIM_INT__s_d_s_s[n * 6 + 3] = SIMINT_FMADD( vrr_const_1_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_d_s_s[n * 6 + 3]);
299
300 PRIM_INT__s_d_s_s[n * 6 + 4] = SIMINT_MUL(P_PB[2], PRIM_INT__s_p_s_s[n * 3 + 1]);
301 PRIM_INT__s_d_s_s[n * 6 + 4] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_d_s_s[n * 6 + 4]);
302
303 PRIM_INT__s_d_s_s[n * 6 + 5] = SIMINT_MUL(P_PB[2], PRIM_INT__s_p_s_s[n * 3 + 2]);
304 PRIM_INT__s_d_s_s[n * 6 + 5] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_d_s_s[n * 6 + 5]);
305 PRIM_INT__s_d_s_s[n * 6 + 5] = SIMINT_FMADD( vrr_const_1_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_d_s_s[n * 6 + 5]);
306
307 }
308
309
310
311 // Forming PRIM_INT__s_f_s_s[7 * 10];
312 for(n = 0; n < 7; ++n) // loop over orders of auxiliary function
313 {
314
315 PRIM_INT__s_f_s_s[n * 10 + 0] = SIMINT_MUL(P_PB[0], PRIM_INT__s_d_s_s[n * 6 + 0]);
316 PRIM_INT__s_f_s_s[n * 10 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 0], PRIM_INT__s_f_s_s[n * 10 + 0]);
317 PRIM_INT__s_f_s_s[n * 10 + 0] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_p_s_s[n * 3 + 0]), PRIM_INT__s_f_s_s[n * 10 + 0]);
318
319 PRIM_INT__s_f_s_s[n * 10 + 1] = SIMINT_MUL(P_PB[1], PRIM_INT__s_d_s_s[n * 6 + 0]);
320 PRIM_INT__s_f_s_s[n * 10 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 0], PRIM_INT__s_f_s_s[n * 10 + 1]);
321
322 PRIM_INT__s_f_s_s[n * 10 + 2] = SIMINT_MUL(P_PB[2], PRIM_INT__s_d_s_s[n * 6 + 0]);
323 PRIM_INT__s_f_s_s[n * 10 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 0], PRIM_INT__s_f_s_s[n * 10 + 2]);
324
325 PRIM_INT__s_f_s_s[n * 10 + 3] = SIMINT_MUL(P_PB[0], PRIM_INT__s_d_s_s[n * 6 + 3]);
326 PRIM_INT__s_f_s_s[n * 10 + 3] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 3], PRIM_INT__s_f_s_s[n * 10 + 3]);
327
328 PRIM_INT__s_f_s_s[n * 10 + 4] = SIMINT_MUL(P_PB[2], PRIM_INT__s_d_s_s[n * 6 + 1]);
329 PRIM_INT__s_f_s_s[n * 10 + 4] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 1], PRIM_INT__s_f_s_s[n * 10 + 4]);
330
331 PRIM_INT__s_f_s_s[n * 10 + 5] = SIMINT_MUL(P_PB[0], PRIM_INT__s_d_s_s[n * 6 + 5]);
332 PRIM_INT__s_f_s_s[n * 10 + 5] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 5], PRIM_INT__s_f_s_s[n * 10 + 5]);
333
334 PRIM_INT__s_f_s_s[n * 10 + 6] = SIMINT_MUL(P_PB[1], PRIM_INT__s_d_s_s[n * 6 + 3]);
335 PRIM_INT__s_f_s_s[n * 10 + 6] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 3], PRIM_INT__s_f_s_s[n * 10 + 6]);
336 PRIM_INT__s_f_s_s[n * 10 + 6] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_p_s_s[n * 3 + 1]), PRIM_INT__s_f_s_s[n * 10 + 6]);
337
338 PRIM_INT__s_f_s_s[n * 10 + 7] = SIMINT_MUL(P_PB[2], PRIM_INT__s_d_s_s[n * 6 + 3]);
339 PRIM_INT__s_f_s_s[n * 10 + 7] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 3], PRIM_INT__s_f_s_s[n * 10 + 7]);
340
341 PRIM_INT__s_f_s_s[n * 10 + 8] = SIMINT_MUL(P_PB[1], PRIM_INT__s_d_s_s[n * 6 + 5]);
342 PRIM_INT__s_f_s_s[n * 10 + 8] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 5], PRIM_INT__s_f_s_s[n * 10 + 8]);
343
344 PRIM_INT__s_f_s_s[n * 10 + 9] = SIMINT_MUL(P_PB[2], PRIM_INT__s_d_s_s[n * 6 + 5]);
345 PRIM_INT__s_f_s_s[n * 10 + 9] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 5], PRIM_INT__s_f_s_s[n * 10 + 9]);
346 PRIM_INT__s_f_s_s[n * 10 + 9] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_p_s_s[n * 3 + 2]), PRIM_INT__s_f_s_s[n * 10 + 9]);
347
348 }
349
350
351 VRR_L_s_f_s_p(
352 PRIM_INT__s_f_s_p,
353 PRIM_INT__s_f_s_s,
354 PRIM_INT__s_d_s_s,
355 Q_PB,
356 aoq_PQ,
357 one_over_2pq,
358 4);
359
360
361
362 // Forming PRIM_INT__s_d_s_p[4 * 18];
363 for(n = 0; n < 4; ++n) // loop over orders of auxiliary function
364 {
365
366 PRIM_INT__s_d_s_p[n * 18 + 0] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_d_s_s[n * 6 + 0]);
367 PRIM_INT__s_d_s_p[n * 18 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 0], PRIM_INT__s_d_s_p[n * 18 + 0]);
368 PRIM_INT__s_d_s_p[n * 18 + 0] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_d_s_p[n * 18 + 0]);
369
370 PRIM_INT__s_d_s_p[n * 18 + 1] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_d_s_s[n * 6 + 0]);
371 PRIM_INT__s_d_s_p[n * 18 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 0], PRIM_INT__s_d_s_p[n * 18 + 1]);
372
373 PRIM_INT__s_d_s_p[n * 18 + 2] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_d_s_s[n * 6 + 0]);
374 PRIM_INT__s_d_s_p[n * 18 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 0], PRIM_INT__s_d_s_p[n * 18 + 2]);
375
376 PRIM_INT__s_d_s_p[n * 18 + 3] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_d_s_s[n * 6 + 1]);
377 PRIM_INT__s_d_s_p[n * 18 + 3] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 1], PRIM_INT__s_d_s_p[n * 18 + 3]);
378 PRIM_INT__s_d_s_p[n * 18 + 3] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_d_s_p[n * 18 + 3]);
379
380 PRIM_INT__s_d_s_p[n * 18 + 4] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_d_s_s[n * 6 + 1]);
381 PRIM_INT__s_d_s_p[n * 18 + 4] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 1], PRIM_INT__s_d_s_p[n * 18 + 4]);
382 PRIM_INT__s_d_s_p[n * 18 + 4] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_d_s_p[n * 18 + 4]);
383
384 PRIM_INT__s_d_s_p[n * 18 + 5] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_d_s_s[n * 6 + 1]);
385 PRIM_INT__s_d_s_p[n * 18 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 1], PRIM_INT__s_d_s_p[n * 18 + 5]);
386
387 PRIM_INT__s_d_s_p[n * 18 + 6] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_d_s_s[n * 6 + 2]);
388 PRIM_INT__s_d_s_p[n * 18 + 6] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 2], PRIM_INT__s_d_s_p[n * 18 + 6]);
389 PRIM_INT__s_d_s_p[n * 18 + 6] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_d_s_p[n * 18 + 6]);
390
391 PRIM_INT__s_d_s_p[n * 18 + 7] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_d_s_s[n * 6 + 2]);
392 PRIM_INT__s_d_s_p[n * 18 + 7] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 2], PRIM_INT__s_d_s_p[n * 18 + 7]);
393
394 PRIM_INT__s_d_s_p[n * 18 + 8] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_d_s_s[n * 6 + 2]);
395 PRIM_INT__s_d_s_p[n * 18 + 8] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 2], PRIM_INT__s_d_s_p[n * 18 + 8]);
396 PRIM_INT__s_d_s_p[n * 18 + 8] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_d_s_p[n * 18 + 8]);
397
398 PRIM_INT__s_d_s_p[n * 18 + 9] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_d_s_s[n * 6 + 3]);
399 PRIM_INT__s_d_s_p[n * 18 + 9] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 3], PRIM_INT__s_d_s_p[n * 18 + 9]);
400
401 PRIM_INT__s_d_s_p[n * 18 + 10] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_d_s_s[n * 6 + 3]);
402 PRIM_INT__s_d_s_p[n * 18 + 10] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 3], PRIM_INT__s_d_s_p[n * 18 + 10]);
403 PRIM_INT__s_d_s_p[n * 18 + 10] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_d_s_p[n * 18 + 10]);
404
405 PRIM_INT__s_d_s_p[n * 18 + 11] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_d_s_s[n * 6 + 3]);
406 PRIM_INT__s_d_s_p[n * 18 + 11] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 3], PRIM_INT__s_d_s_p[n * 18 + 11]);
407
408 PRIM_INT__s_d_s_p[n * 18 + 12] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_d_s_s[n * 6 + 4]);
409 PRIM_INT__s_d_s_p[n * 18 + 12] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 4], PRIM_INT__s_d_s_p[n * 18 + 12]);
410
411 PRIM_INT__s_d_s_p[n * 18 + 13] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_d_s_s[n * 6 + 4]);
412 PRIM_INT__s_d_s_p[n * 18 + 13] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 4], PRIM_INT__s_d_s_p[n * 18 + 13]);
413 PRIM_INT__s_d_s_p[n * 18 + 13] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_d_s_p[n * 18 + 13]);
414
415 PRIM_INT__s_d_s_p[n * 18 + 14] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_d_s_s[n * 6 + 4]);
416 PRIM_INT__s_d_s_p[n * 18 + 14] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 4], PRIM_INT__s_d_s_p[n * 18 + 14]);
417 PRIM_INT__s_d_s_p[n * 18 + 14] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_d_s_p[n * 18 + 14]);
418
419 PRIM_INT__s_d_s_p[n * 18 + 15] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_d_s_s[n * 6 + 5]);
420 PRIM_INT__s_d_s_p[n * 18 + 15] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_d_s_s[(n+1) * 6 + 5], PRIM_INT__s_d_s_p[n * 18 + 15]);
421
422 PRIM_INT__s_d_s_p[n * 18 + 16] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_d_s_s[n * 6 + 5]);
423 PRIM_INT__s_d_s_p[n * 18 + 16] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_d_s_s[(n+1) * 6 + 5], PRIM_INT__s_d_s_p[n * 18 + 16]);
424
425 PRIM_INT__s_d_s_p[n * 18 + 17] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_d_s_s[n * 6 + 5]);
426 PRIM_INT__s_d_s_p[n * 18 + 17] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_d_s_s[(n+1) * 6 + 5], PRIM_INT__s_d_s_p[n * 18 + 17]);
427 PRIM_INT__s_d_s_p[n * 18 + 17] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_d_s_p[n * 18 + 17]);
428
429 }
430
431
432 VRR_L_s_f_s_d(
433 PRIM_INT__s_f_s_d,
434 PRIM_INT__s_f_s_p,
435 PRIM_INT__s_f_s_s,
436 PRIM_INT__s_d_s_p,
437 Q_PB,
438 a_over_q,
439 aoq_PQ,
440 one_over_2pq,
441 one_over_2q,
442 3);
443
444
445
446 // Forming PRIM_INT__s_p_s_p[4 * 9];
447 for(n = 0; n < 4; ++n) // loop over orders of auxiliary function
448 {
449
450 PRIM_INT__s_p_s_p[n * 9 + 0] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_p_s_s[n * 3 + 0]);
451 PRIM_INT__s_p_s_p[n * 9 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_p_s_p[n * 9 + 0]);
452 PRIM_INT__s_p_s_p[n * 9 + 0] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_p_s_p[n * 9 + 0]);
453
454 PRIM_INT__s_p_s_p[n * 9 + 1] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_p_s_s[n * 3 + 0]);
455 PRIM_INT__s_p_s_p[n * 9 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_p_s_p[n * 9 + 1]);
456
457 PRIM_INT__s_p_s_p[n * 9 + 2] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_p_s_s[n * 3 + 0]);
458 PRIM_INT__s_p_s_p[n * 9 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_p_s_p[n * 9 + 2]);
459
460 PRIM_INT__s_p_s_p[n * 9 + 3] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_p_s_s[n * 3 + 1]);
461 PRIM_INT__s_p_s_p[n * 9 + 3] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_p_s_p[n * 9 + 3]);
462
463 PRIM_INT__s_p_s_p[n * 9 + 4] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_p_s_s[n * 3 + 1]);
464 PRIM_INT__s_p_s_p[n * 9 + 4] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_p_s_p[n * 9 + 4]);
465 PRIM_INT__s_p_s_p[n * 9 + 4] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_p_s_p[n * 9 + 4]);
466
467 PRIM_INT__s_p_s_p[n * 9 + 5] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_p_s_s[n * 3 + 1]);
468 PRIM_INT__s_p_s_p[n * 9 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_p_s_p[n * 9 + 5]);
469
470 PRIM_INT__s_p_s_p[n * 9 + 6] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_p_s_s[n * 3 + 2]);
471 PRIM_INT__s_p_s_p[n * 9 + 6] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_p_s_p[n * 9 + 6]);
472
473 PRIM_INT__s_p_s_p[n * 9 + 7] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_p_s_s[n * 3 + 2]);
474 PRIM_INT__s_p_s_p[n * 9 + 7] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_p_s_p[n * 9 + 7]);
475
476 PRIM_INT__s_p_s_p[n * 9 + 8] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_p_s_s[n * 3 + 2]);
477 PRIM_INT__s_p_s_p[n * 9 + 8] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_p_s_p[n * 9 + 8]);
478 PRIM_INT__s_p_s_p[n * 9 + 8] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_p_s_p[n * 9 + 8]);
479
480 }
481
482
483 VRR_L_s_d_s_d(
484 PRIM_INT__s_d_s_d,
485 PRIM_INT__s_d_s_p,
486 PRIM_INT__s_d_s_s,
487 PRIM_INT__s_p_s_p,
488 Q_PB,
489 a_over_q,
490 aoq_PQ,
491 one_over_2pq,
492 one_over_2q,
493 3);
494
495
496 ostei_general_vrr_L(0, 3, 0, 3, 2,
497 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
498 PRIM_INT__s_f_s_d, NULL, PRIM_INT__s_f_s_p, NULL, PRIM_INT__s_d_s_d, PRIM_INT__s_f_s_f);
499
500
501
502 // Forming PRIM_INT__s_s_s_p[4 * 3];
503 for(n = 0; n < 4; ++n) // loop over orders of auxiliary function
504 {
505
506 PRIM_INT__s_s_s_p[n * 3 + 0] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
507 PRIM_INT__s_s_s_p[n * 3 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_p[n * 3 + 0]);
508
509 PRIM_INT__s_s_s_p[n * 3 + 1] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
510 PRIM_INT__s_s_s_p[n * 3 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_p[n * 3 + 1]);
511
512 PRIM_INT__s_s_s_p[n * 3 + 2] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
513 PRIM_INT__s_s_s_p[n * 3 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_p[n * 3 + 2]);
514
515 }
516
517
518
519 // Forming PRIM_INT__s_p_s_d[3 * 18];
520 for(n = 0; n < 3; ++n) // loop over orders of auxiliary function
521 {
522
523 PRIM_INT__s_p_s_d[n * 18 + 0] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_p_s_p[n * 9 + 0]);
524 PRIM_INT__s_p_s_d[n * 18 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_p_s_p[(n+1) * 9 + 0], PRIM_INT__s_p_s_d[n * 18 + 0]);
525 PRIM_INT__s_p_s_d[n * 18 + 0] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_p_s_s[n * 3 + 0]), PRIM_INT__s_p_s_d[n * 18 + 0]);
526 PRIM_INT__s_p_s_d[n * 18 + 0] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_s_s_p[(n+1) * 3 + 0], PRIM_INT__s_p_s_d[n * 18 + 0]);
527
528 PRIM_INT__s_p_s_d[n * 18 + 3] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_p_s_p[n * 9 + 1]);
529 PRIM_INT__s_p_s_d[n * 18 + 3] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_p_s_p[(n+1) * 9 + 1], PRIM_INT__s_p_s_d[n * 18 + 3]);
530 PRIM_INT__s_p_s_d[n * 18 + 3] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_p_s_s[n * 3 + 0]), PRIM_INT__s_p_s_d[n * 18 + 3]);
531
532 PRIM_INT__s_p_s_d[n * 18 + 5] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_p_s_p[n * 9 + 2]);
533 PRIM_INT__s_p_s_d[n * 18 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_p_s_p[(n+1) * 9 + 2], PRIM_INT__s_p_s_d[n * 18 + 5]);
534 PRIM_INT__s_p_s_d[n * 18 + 5] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 0], PRIM_INT__s_p_s_s[n * 3 + 0]), PRIM_INT__s_p_s_d[n * 18 + 5]);
535
536 PRIM_INT__s_p_s_d[n * 18 + 6] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_p_s_p[n * 9 + 3]);
537 PRIM_INT__s_p_s_d[n * 18 + 6] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_p_s_p[(n+1) * 9 + 3], PRIM_INT__s_p_s_d[n * 18 + 6]);
538 PRIM_INT__s_p_s_d[n * 18 + 6] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_p_s_s[n * 3 + 1]), PRIM_INT__s_p_s_d[n * 18 + 6]);
539
540 PRIM_INT__s_p_s_d[n * 18 + 9] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_p_s_p[n * 9 + 4]);
541 PRIM_INT__s_p_s_d[n * 18 + 9] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_p_s_p[(n+1) * 9 + 4], PRIM_INT__s_p_s_d[n * 18 + 9]);
542 PRIM_INT__s_p_s_d[n * 18 + 9] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_p_s_s[n * 3 + 1]), PRIM_INT__s_p_s_d[n * 18 + 9]);
543 PRIM_INT__s_p_s_d[n * 18 + 9] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_s_s_p[(n+1) * 3 + 1], PRIM_INT__s_p_s_d[n * 18 + 9]);
544
545 PRIM_INT__s_p_s_d[n * 18 + 11] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_p_s_p[n * 9 + 5]);
546 PRIM_INT__s_p_s_d[n * 18 + 11] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_p_s_p[(n+1) * 9 + 5], PRIM_INT__s_p_s_d[n * 18 + 11]);
547 PRIM_INT__s_p_s_d[n * 18 + 11] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 1], PRIM_INT__s_p_s_s[n * 3 + 1]), PRIM_INT__s_p_s_d[n * 18 + 11]);
548
549 PRIM_INT__s_p_s_d[n * 18 + 12] = SIMINT_MUL(Q_PB[0], PRIM_INT__s_p_s_p[n * 9 + 6]);
550 PRIM_INT__s_p_s_d[n * 18 + 12] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_p_s_p[(n+1) * 9 + 6], PRIM_INT__s_p_s_d[n * 18 + 12]);
551 PRIM_INT__s_p_s_d[n * 18 + 12] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_p_s_s[n * 3 + 2]), PRIM_INT__s_p_s_d[n * 18 + 12]);
552
553 PRIM_INT__s_p_s_d[n * 18 + 15] = SIMINT_MUL(Q_PB[1], PRIM_INT__s_p_s_p[n * 9 + 7]);
554 PRIM_INT__s_p_s_d[n * 18 + 15] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_p_s_p[(n+1) * 9 + 7], PRIM_INT__s_p_s_d[n * 18 + 15]);
555 PRIM_INT__s_p_s_d[n * 18 + 15] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_p_s_s[n * 3 + 2]), PRIM_INT__s_p_s_d[n * 18 + 15]);
556
557 PRIM_INT__s_p_s_d[n * 18 + 17] = SIMINT_MUL(Q_PB[2], PRIM_INT__s_p_s_p[n * 9 + 8]);
558 PRIM_INT__s_p_s_d[n * 18 + 17] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_p_s_p[(n+1) * 9 + 8], PRIM_INT__s_p_s_d[n * 18 + 17]);
559 PRIM_INT__s_p_s_d[n * 18 + 17] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_p_s_s[(n+1) * 3 + 2], PRIM_INT__s_p_s_s[n * 3 + 2]), PRIM_INT__s_p_s_d[n * 18 + 17]);
560 PRIM_INT__s_p_s_d[n * 18 + 17] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__s_s_s_p[(n+1) * 3 + 2], PRIM_INT__s_p_s_d[n * 18 + 17]);
561
562 }
563
564
565 VRR_L_s_d_s_f(
566 PRIM_INT__s_d_s_f,
567 PRIM_INT__s_d_s_d,
568 PRIM_INT__s_d_s_p,
569 PRIM_INT__s_p_s_d,
570 Q_PB,
571 a_over_q,
572 aoq_PQ,
573 one_over_2pq,
574 one_over_2q,
575 2);
576
577
578 ostei_general_vrr_L(0, 3, 0, 4, 1,
579 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
580 PRIM_INT__s_f_s_f, NULL, PRIM_INT__s_f_s_d, NULL, PRIM_INT__s_d_s_f, PRIM_INT__s_f_s_g);
581
582
583 VRR_J_s_g_s_s(
584 PRIM_INT__s_g_s_s,
585 PRIM_INT__s_f_s_s,
586 PRIM_INT__s_d_s_s,
587 P_PB,
588 a_over_p,
589 aop_PQ,
590 one_over_2p,
591 6);
592
593
594 VRR_L_s_g_s_p(
595 PRIM_INT__s_g_s_p,
596 PRIM_INT__s_g_s_s,
597 PRIM_INT__s_f_s_s,
598 Q_PB,
599 aoq_PQ,
600 one_over_2pq,
601 4);
602
603
604 ostei_general_vrr_L(0, 4, 0, 2, 3,
605 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
606 PRIM_INT__s_g_s_p, NULL, PRIM_INT__s_g_s_s, NULL, PRIM_INT__s_f_s_p, PRIM_INT__s_g_s_d);
607
608
609 ostei_general_vrr_L(0, 4, 0, 3, 2,
610 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
611 PRIM_INT__s_g_s_d, NULL, PRIM_INT__s_g_s_p, NULL, PRIM_INT__s_f_s_d, PRIM_INT__s_g_s_f);
612
613
614 ostei_general_vrr_L(0, 4, 0, 4, 1,
615 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
616 PRIM_INT__s_g_s_f, NULL, PRIM_INT__s_g_s_d, NULL, PRIM_INT__s_f_s_f, PRIM_INT__s_g_s_g);
617
618
619 VRR_J_s_h_s_s(
620 PRIM_INT__s_h_s_s,
621 PRIM_INT__s_g_s_s,
622 PRIM_INT__s_f_s_s,
623 P_PB,
624 a_over_p,
625 aop_PQ,
626 one_over_2p,
627 5);
628
629
630 ostei_general_vrr_L(0, 5, 0, 1, 4,
631 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
632 PRIM_INT__s_h_s_s, NULL, NULL, NULL, PRIM_INT__s_g_s_s, PRIM_INT__s_h_s_p);
633
634
635 ostei_general_vrr_L(0, 5, 0, 2, 3,
636 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
637 PRIM_INT__s_h_s_p, NULL, PRIM_INT__s_h_s_s, NULL, PRIM_INT__s_g_s_p, PRIM_INT__s_h_s_d);
638
639
640 ostei_general_vrr_L(0, 5, 0, 3, 2,
641 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
642 PRIM_INT__s_h_s_d, NULL, PRIM_INT__s_h_s_p, NULL, PRIM_INT__s_g_s_d, PRIM_INT__s_h_s_f);
643
644
645 ostei_general_vrr_L(0, 5, 0, 4, 1,
646 one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PB,
647 PRIM_INT__s_h_s_f, NULL, PRIM_INT__s_h_s_d, NULL, PRIM_INT__s_g_s_f, PRIM_INT__s_h_s_g);
648
649
650
651
652 ////////////////////////////////////
653 // Accumulate contracted integrals
654 ////////////////////////////////////
655 if(lastoffset == 0)
656 {
657 contract_all(100, PRIM_INT__s_f_s_f, PRIM_PTR_INT__s_f_s_f);
658 contract_all(150, PRIM_INT__s_f_s_g, PRIM_PTR_INT__s_f_s_g);
659 contract_all(150, PRIM_INT__s_g_s_f, PRIM_PTR_INT__s_g_s_f);
660 contract_all(225, PRIM_INT__s_g_s_g, PRIM_PTR_INT__s_g_s_g);
661 contract_all(210, PRIM_INT__s_h_s_f, PRIM_PTR_INT__s_h_s_f);
662 contract_all(315, PRIM_INT__s_h_s_g, PRIM_PTR_INT__s_h_s_g);
663 }
664 else
665 {
666 contract(100, shelloffsets, PRIM_INT__s_f_s_f, PRIM_PTR_INT__s_f_s_f);
667 contract(150, shelloffsets, PRIM_INT__s_f_s_g, PRIM_PTR_INT__s_f_s_g);
668 contract(150, shelloffsets, PRIM_INT__s_g_s_f, PRIM_PTR_INT__s_g_s_f);
669 contract(225, shelloffsets, PRIM_INT__s_g_s_g, PRIM_PTR_INT__s_g_s_g);
670 contract(210, shelloffsets, PRIM_INT__s_h_s_f, PRIM_PTR_INT__s_h_s_f);
671 contract(315, shelloffsets, PRIM_INT__s_h_s_g, PRIM_PTR_INT__s_h_s_g);
672 PRIM_PTR_INT__s_f_s_f += lastoffset*100;
673 PRIM_PTR_INT__s_f_s_g += lastoffset*150;
674 PRIM_PTR_INT__s_g_s_f += lastoffset*150;
675 PRIM_PTR_INT__s_g_s_g += lastoffset*225;
676 PRIM_PTR_INT__s_h_s_f += lastoffset*210;
677 PRIM_PTR_INT__s_h_s_g += lastoffset*315;
678 }
679
680 } // close loop over j
681 } // close loop over i
682
683 //Advance to the next batch
684 jstart = SIMINT_SIMD_ROUND(jend);
685
686 //////////////////////////////////////////////
687 // Contracted integrals: Horizontal recurrance
688 //////////////////////////////////////////////
689
690
691 const double hAB[3] = { P.AB_x[ab], P.AB_y[ab], P.AB_z[ab] };
692
693
694 for(abcd = 0; abcd < nshellbatch; ++abcd, ++real_abcd)
695 {
696 const double hCD[3] = { Q.AB_x[cd+abcd], Q.AB_y[cd+abcd], Q.AB_z[cd+abcd] };
697
698 // set up HRR pointers
699 double const * restrict HRR_INT__s_f_s_f = INT__s_f_s_f + abcd * 100;
700 double const * restrict HRR_INT__s_f_s_g = INT__s_f_s_g + abcd * 150;
701 double const * restrict HRR_INT__s_g_s_f = INT__s_g_s_f + abcd * 150;
702 double const * restrict HRR_INT__s_g_s_g = INT__s_g_s_g + abcd * 225;
703 double const * restrict HRR_INT__s_h_s_f = INT__s_h_s_f + abcd * 210;
704 double const * restrict HRR_INT__s_h_s_g = INT__s_h_s_g + abcd * 315;
705 double * restrict HRR_INT__d_f_p_f = INT__d_f_p_f + real_abcd * 1800;
706
707 // form INT__p_f_s_f
708 HRR_I_p_f(
709 HRR_INT__p_f_s_f,
710 HRR_INT__s_f_s_f,
711 HRR_INT__s_g_s_f,
712 hAB, 10);
713
714 // form INT__p_f_s_g
715 HRR_I_p_f(
716 HRR_INT__p_f_s_g,
717 HRR_INT__s_f_s_g,
718 HRR_INT__s_g_s_g,
719 hAB, 15);
720
721 // form INT__p_g_s_f
722 HRR_I_p_g(
723 HRR_INT__p_g_s_f,
724 HRR_INT__s_g_s_f,
725 HRR_INT__s_h_s_f,
726 hAB, 10);
727
728 // form INT__p_g_s_g
729 HRR_I_p_g(
730 HRR_INT__p_g_s_g,
731 HRR_INT__s_g_s_g,
732 HRR_INT__s_h_s_g,
733 hAB, 15);
734
735 // form INT__d_f_s_f
736 HRR_I_d_f(
737 HRR_INT__d_f_s_f,
738 HRR_INT__p_f_s_f,
739 HRR_INT__p_g_s_f,
740 hAB, 10);
741
742 // form INT__d_f_s_g
743 HRR_I_d_f(
744 HRR_INT__d_f_s_g,
745 HRR_INT__p_f_s_g,
746 HRR_INT__p_g_s_g,
747 hAB, 15);
748
749 // form INT__d_f_p_f
750 HRR_K_p_f(
751 HRR_INT__d_f_p_f,
752 HRR_INT__d_f_s_f,
753 HRR_INT__d_f_s_g,
754 hCD, 60);
755
756
757 } // close HRR loop
758
759
760 } // close loop cdbatch
761
762 istart = iend;
763 } // close loop over ab
764
765 return P.nshell12_clip * Q.nshell12_clip;
766 }
767
768