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