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