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_h_h_h_p(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__h_h_h_p)8 int ostei_h_h_h_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__h_h_h_p)
13 {
14 
15     SIMINT_ASSUME_ALIGN_DBL(work);
16     SIMINT_ASSUME_ALIGN_DBL(INT__h_h_h_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__h_s_h_s = work + (SIMINT_NSHELL_SIMD * 0);
30     double * const INT__h_s_i_s = work + (SIMINT_NSHELL_SIMD * 441);
31     double * const INT__i_s_h_s = work + (SIMINT_NSHELL_SIMD * 1029);
32     double * const INT__i_s_i_s = work + (SIMINT_NSHELL_SIMD * 1617);
33     double * const INT__k_s_h_s = work + (SIMINT_NSHELL_SIMD * 2401);
34     double * const INT__k_s_i_s = work + (SIMINT_NSHELL_SIMD * 3157);
35     double * const INT__l_s_h_s = work + (SIMINT_NSHELL_SIMD * 4165);
36     double * const INT__l_s_i_s = work + (SIMINT_NSHELL_SIMD * 5110);
37     double * const INT__m_s_h_s = work + (SIMINT_NSHELL_SIMD * 6370);
38     double * const INT__m_s_i_s = work + (SIMINT_NSHELL_SIMD * 7525);
39     double * const INT__n_s_h_s = work + (SIMINT_NSHELL_SIMD * 9065);
40     double * const INT__n_s_i_s = work + (SIMINT_NSHELL_SIMD * 10451);
41     SIMINT_DBLTYPE * const primwork = (SIMINT_DBLTYPE *)(work + SIMINT_NSHELL_SIMD*12299);
42     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_s = primwork + 0;
43     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_p_s = primwork + 17;
44     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_s_s = primwork + 35;
45     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_p_s = primwork + 83;
46     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_d_s = primwork + 137;
47     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_s_s = primwork + 227;
48     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_p_s = primwork + 317;
49     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_d_s = primwork + 425;
50     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_f_s = primwork + 605;
51     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_s_s = primwork + 845;
52     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_p_s = primwork + 985;
53     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_d_s = primwork + 1165;
54     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_f_s = primwork + 1465;
55     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_g_s = primwork + 1865;
56     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_s_s = primwork + 2315;
57     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_p_s = primwork + 2510;
58     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_d_s = primwork + 2780;
59     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_f_s = primwork + 3230;
60     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_g_s = primwork + 3830;
61     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_h_s = primwork + 4505;
62     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_s_s = primwork + 5135;
63     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_p_s = primwork + 5387;
64     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_d_s = primwork + 5765;
65     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_f_s = primwork + 6395;
66     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_g_s = primwork + 7235;
67     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_h_s = primwork + 8180;
68     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_i_s = primwork + 9062;
69     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_s_s = primwork + 9650;
70     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_p_s = primwork + 9958;
71     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_d_s = primwork + 10462;
72     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_f_s = primwork + 11302;
73     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_g_s = primwork + 12422;
74     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_h_s = primwork + 13682;
75     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_i_s = primwork + 14858;
76     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_s_s = primwork + 15642;
77     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_p_s = primwork + 16002;
78     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_d_s = primwork + 16650;
79     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_f_s = primwork + 17730;
80     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_g_s = primwork + 19170;
81     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_h_s = primwork + 20790;
82     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_i_s = primwork + 22302;
83     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_s_s = primwork + 23310;
84     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_p_s = primwork + 23715;
85     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_d_s = primwork + 24525;
86     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_f_s = primwork + 25875;
87     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_g_s = primwork + 27675;
88     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_h_s = primwork + 29700;
89     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_i_s = primwork + 31590;
90     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_s_s = primwork + 32850;
91     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_p_s = primwork + 33290;
92     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_d_s = primwork + 34280;
93     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_f_s = primwork + 35930;
94     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_g_s = primwork + 38130;
95     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_h_s = primwork + 40605;
96     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_i_s = primwork + 42915;
97     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_s_s = primwork + 44455;
98     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_p_s = primwork + 44917;
99     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_d_s = primwork + 46105;
100     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_f_s = primwork + 48085;
101     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_g_s = primwork + 50725;
102     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_h_s = primwork + 53695;
103     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_i_s = primwork + 56467;
104     double * const hrrwork = (double *)(primwork + 58315);
105     double * const HRR_INT__h_p_h_s = hrrwork + 0;
106     double * const HRR_INT__h_p_i_s = hrrwork + 1323;
107     double * const HRR_INT__h_d_h_s = hrrwork + 3087;
108     double * const HRR_INT__h_d_i_s = hrrwork + 5733;
109     double * const HRR_INT__h_f_h_s = hrrwork + 9261;
110     double * const HRR_INT__h_f_i_s = hrrwork + 13671;
111     double * const HRR_INT__h_g_h_s = hrrwork + 19551;
112     double * const HRR_INT__h_g_i_s = hrrwork + 26166;
113     double * const HRR_INT__h_h_h_s = hrrwork + 34986;
114     double * const HRR_INT__h_h_i_s = hrrwork + 44247;
115     double * const HRR_INT__i_p_h_s = hrrwork + 56595;
116     double * const HRR_INT__i_p_i_s = hrrwork + 58359;
117     double * const HRR_INT__i_d_h_s = hrrwork + 60711;
118     double * const HRR_INT__i_d_i_s = hrrwork + 64239;
119     double * const HRR_INT__i_f_h_s = hrrwork + 68943;
120     double * const HRR_INT__i_f_i_s = hrrwork + 74823;
121     double * const HRR_INT__i_g_h_s = hrrwork + 82663;
122     double * const HRR_INT__i_g_i_s = hrrwork + 91483;
123     double * const HRR_INT__k_p_h_s = hrrwork + 103243;
124     double * const HRR_INT__k_p_i_s = hrrwork + 105511;
125     double * const HRR_INT__k_d_h_s = hrrwork + 108535;
126     double * const HRR_INT__k_d_i_s = hrrwork + 113071;
127     double * const HRR_INT__k_f_h_s = hrrwork + 119119;
128     double * const HRR_INT__k_f_i_s = hrrwork + 126679;
129     double * const HRR_INT__l_p_h_s = hrrwork + 136759;
130     double * const HRR_INT__l_p_i_s = hrrwork + 139594;
131     double * const HRR_INT__l_d_h_s = hrrwork + 143374;
132     double * const HRR_INT__l_d_i_s = hrrwork + 149044;
133     double * const HRR_INT__m_p_h_s = hrrwork + 156604;
134     double * const HRR_INT__m_p_i_s = hrrwork + 160069;
135 
136 
137     // Create constants
138     const SIMINT_DBLTYPE const_1 = SIMINT_DBLSET1(1);
139     const SIMINT_DBLTYPE const_10 = SIMINT_DBLSET1(10);
140     const SIMINT_DBLTYPE const_2 = SIMINT_DBLSET1(2);
141     const SIMINT_DBLTYPE const_3 = SIMINT_DBLSET1(3);
142     const SIMINT_DBLTYPE const_4 = SIMINT_DBLSET1(4);
143     const SIMINT_DBLTYPE const_5 = SIMINT_DBLSET1(5);
144     const SIMINT_DBLTYPE const_6 = SIMINT_DBLSET1(6);
145     const SIMINT_DBLTYPE const_7 = SIMINT_DBLSET1(7);
146     const SIMINT_DBLTYPE const_8 = SIMINT_DBLSET1(8);
147     const SIMINT_DBLTYPE const_9 = SIMINT_DBLSET1(9);
148     const SIMINT_DBLTYPE one_half = SIMINT_DBLSET1(0.5);
149 
150 
151     ////////////////////////////////////////
152     // Loop over shells and primitives
153     ////////////////////////////////////////
154 
155     real_abcd = 0;
156     istart = 0;
157     for(ab = 0; ab < P.nshell12_clip; ++ab)
158     {
159         const int iend = istart + P.nprim12[ab];
160 
161         cd = 0;
162         jstart = 0;
163 
164         for(cd = 0; cd < Q.nshell12_clip; cd += SIMINT_NSHELL_SIMD)
165         {
166             const int nshellbatch = ((cd + SIMINT_NSHELL_SIMD) > Q.nshell12_clip) ? Q.nshell12_clip - cd : SIMINT_NSHELL_SIMD;
167             int jend = jstart;
168             for(i = 0; i < nshellbatch; i++)
169                 jend += Q.nprim12[cd+i];
170 
171             // Clear the beginning of the workspace (where we are accumulating integrals)
172             memset(work, 0, SIMINT_NSHELL_SIMD * 12299 * sizeof(double));
173             abcd = 0;
174 
175 
176             for(i = istart; i < iend; ++i)
177             {
178                 SIMINT_DBLTYPE bra_screen_max;  // only used if check_screen
179 
180                 if(check_screen)
181                 {
182                     // Skip this whole thing if always insignificant
183                     if((P.screen[i] * Q.screen_max) < screen_tol)
184                         continue;
185                     bra_screen_max = SIMINT_DBLSET1(P.screen[i]);
186                 }
187 
188                 icd = 0;
189                 iprimcd = 0;
190                 nprim_icd = Q.nprim12[cd];
191                 double * restrict PRIM_PTR_INT__h_s_h_s = INT__h_s_h_s + abcd * 441;
192                 double * restrict PRIM_PTR_INT__h_s_i_s = INT__h_s_i_s + abcd * 588;
193                 double * restrict PRIM_PTR_INT__i_s_h_s = INT__i_s_h_s + abcd * 588;
194                 double * restrict PRIM_PTR_INT__i_s_i_s = INT__i_s_i_s + abcd * 784;
195                 double * restrict PRIM_PTR_INT__k_s_h_s = INT__k_s_h_s + abcd * 756;
196                 double * restrict PRIM_PTR_INT__k_s_i_s = INT__k_s_i_s + abcd * 1008;
197                 double * restrict PRIM_PTR_INT__l_s_h_s = INT__l_s_h_s + abcd * 945;
198                 double * restrict PRIM_PTR_INT__l_s_i_s = INT__l_s_i_s + abcd * 1260;
199                 double * restrict PRIM_PTR_INT__m_s_h_s = INT__m_s_h_s + abcd * 1155;
200                 double * restrict PRIM_PTR_INT__m_s_i_s = INT__m_s_i_s + abcd * 1540;
201                 double * restrict PRIM_PTR_INT__n_s_h_s = INT__n_s_h_s + abcd * 1386;
202                 double * restrict PRIM_PTR_INT__n_s_i_s = INT__n_s_i_s + abcd * 1848;
203 
204 
205 
206                 // Load these one per loop over i
207                 const SIMINT_DBLTYPE P_alpha = SIMINT_DBLSET1(P.alpha[i]);
208                 const SIMINT_DBLTYPE P_prefac = SIMINT_DBLSET1(P.prefac[i]);
209                 const SIMINT_DBLTYPE Pxyz[3] = { SIMINT_DBLSET1(P.x[i]), SIMINT_DBLSET1(P.y[i]), SIMINT_DBLSET1(P.z[i]) };
210 
211                 const SIMINT_DBLTYPE P_PA[3] = { SIMINT_DBLSET1(P.PA_x[i]), SIMINT_DBLSET1(P.PA_y[i]), SIMINT_DBLSET1(P.PA_z[i]) };
212 
213                 for(j = jstart; j < jend; j += SIMINT_SIMD_LEN)
214                 {
215                     // calculate the shell offsets
216                     // these are the offset from the shell pointed to by cd
217                     // for each element
218                     int shelloffsets[SIMINT_SIMD_LEN] = {0};
219                     int lastoffset = 0;
220                     const int nlane = ( ((j + SIMINT_SIMD_LEN) < jend) ? SIMINT_SIMD_LEN : (jend - j));
221 
222                     if((iprimcd + SIMINT_SIMD_LEN) >= nprim_icd)
223                     {
224                         // Handle if the first element of the vector is a new shell
225                         if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
226                         {
227                             nprim_icd += Q.nprim12[cd + (++icd)];
228                             PRIM_PTR_INT__h_s_h_s += 441;
229                             PRIM_PTR_INT__h_s_i_s += 588;
230                             PRIM_PTR_INT__i_s_h_s += 588;
231                             PRIM_PTR_INT__i_s_i_s += 784;
232                             PRIM_PTR_INT__k_s_h_s += 756;
233                             PRIM_PTR_INT__k_s_i_s += 1008;
234                             PRIM_PTR_INT__l_s_h_s += 945;
235                             PRIM_PTR_INT__l_s_i_s += 1260;
236                             PRIM_PTR_INT__m_s_h_s += 1155;
237                             PRIM_PTR_INT__m_s_i_s += 1540;
238                             PRIM_PTR_INT__n_s_h_s += 1386;
239                             PRIM_PTR_INT__n_s_i_s += 1848;
240                         }
241                         iprimcd++;
242                         for(n = 1; n < SIMINT_SIMD_LEN; ++n)
243                         {
244                             if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
245                             {
246                                 shelloffsets[n] = shelloffsets[n-1] + 1;
247                                 lastoffset++;
248                                 nprim_icd += Q.nprim12[cd + (++icd)];
249                             }
250                             else
251                                 shelloffsets[n] = shelloffsets[n-1];
252                             iprimcd++;
253                         }
254                     }
255                     else
256                         iprimcd += SIMINT_SIMD_LEN;
257 
258                     // Do we have to compute this vector (or has it been screened out)?
259                     // (not_screened != 0 means we have to do this vector)
260                     if(check_screen)
261                     {
262                         const double vmax = vector_max(SIMINT_MUL(bra_screen_max, SIMINT_DBLLOAD(Q.screen, j)));
263                         if(vmax < screen_tol)
264                         {
265                             PRIM_PTR_INT__h_s_h_s += lastoffset*441;
266                             PRIM_PTR_INT__h_s_i_s += lastoffset*588;
267                             PRIM_PTR_INT__i_s_h_s += lastoffset*588;
268                             PRIM_PTR_INT__i_s_i_s += lastoffset*784;
269                             PRIM_PTR_INT__k_s_h_s += lastoffset*756;
270                             PRIM_PTR_INT__k_s_i_s += lastoffset*1008;
271                             PRIM_PTR_INT__l_s_h_s += lastoffset*945;
272                             PRIM_PTR_INT__l_s_i_s += lastoffset*1260;
273                             PRIM_PTR_INT__m_s_h_s += lastoffset*1155;
274                             PRIM_PTR_INT__m_s_i_s += lastoffset*1540;
275                             PRIM_PTR_INT__n_s_h_s += lastoffset*1386;
276                             PRIM_PTR_INT__n_s_i_s += lastoffset*1848;
277                             continue;
278                         }
279                     }
280 
281                     const SIMINT_DBLTYPE Q_alpha = SIMINT_DBLLOAD(Q.alpha, j);
282                     const SIMINT_DBLTYPE PQalpha_mul = SIMINT_MUL(P_alpha, Q_alpha);
283                     const SIMINT_DBLTYPE PQalpha_sum = SIMINT_ADD(P_alpha, Q_alpha);
284                     const SIMINT_DBLTYPE one_over_PQalpha_sum = SIMINT_DIV(const_1, PQalpha_sum);
285 
286 
287                     /* construct R2 = (Px - Qx)**2 + (Py - Qy)**2 + (Pz -Qz)**2 */
288                     SIMINT_DBLTYPE PQ[3];
289                     PQ[0] = SIMINT_SUB(Pxyz[0], SIMINT_DBLLOAD(Q.x, j));
290                     PQ[1] = SIMINT_SUB(Pxyz[1], SIMINT_DBLLOAD(Q.y, j));
291                     PQ[2] = SIMINT_SUB(Pxyz[2], SIMINT_DBLLOAD(Q.z, j));
292                     SIMINT_DBLTYPE R2 = SIMINT_MUL(PQ[0], PQ[0]);
293                     R2 = SIMINT_FMADD(PQ[1], PQ[1], R2);
294                     R2 = SIMINT_FMADD(PQ[2], PQ[2], R2);
295 
296                     const SIMINT_DBLTYPE alpha = SIMINT_MUL(PQalpha_mul, one_over_PQalpha_sum); // alpha from MEST
297                     const SIMINT_DBLTYPE one_over_p = SIMINT_DIV(const_1, P_alpha);
298                     const SIMINT_DBLTYPE one_over_q = SIMINT_DIV(const_1, Q_alpha);
299                     const SIMINT_DBLTYPE one_over_2p = SIMINT_MUL(one_half, one_over_p);
300                     const SIMINT_DBLTYPE one_over_2q = SIMINT_MUL(one_half, one_over_q);
301                     const SIMINT_DBLTYPE one_over_2pq = SIMINT_MUL(one_half, one_over_PQalpha_sum);
302                     const SIMINT_DBLTYPE Q_PA[3] = { SIMINT_DBLLOAD(Q.PA_x, j), SIMINT_DBLLOAD(Q.PA_y, j), SIMINT_DBLLOAD(Q.PA_z, j) };
303 
304                     // NOTE: Minus sign!
305                     const SIMINT_DBLTYPE a_over_p = SIMINT_MUL(SIMINT_NEG(alpha), one_over_p);
306                     SIMINT_DBLTYPE aop_PQ[3];
307                     aop_PQ[0] = SIMINT_MUL(a_over_p, PQ[0]);
308                     aop_PQ[1] = SIMINT_MUL(a_over_p, PQ[1]);
309                     aop_PQ[2] = SIMINT_MUL(a_over_p, PQ[2]);
310 
311                     SIMINT_DBLTYPE a_over_q = SIMINT_MUL(alpha, one_over_q);
312                     SIMINT_DBLTYPE aoq_PQ[3];
313                     aoq_PQ[0] = SIMINT_MUL(a_over_q, PQ[0]);
314                     aoq_PQ[1] = SIMINT_MUL(a_over_q, PQ[1]);
315                     aoq_PQ[2] = SIMINT_MUL(a_over_q, PQ[2]);
316                     // Put a minus sign here so we don't have to in RR routines
317                     a_over_q = SIMINT_NEG(a_over_q);
318 
319 
320                     //////////////////////////////////////////////
321                     // Fjt function section
322                     // Maximum v value: 16
323                     //////////////////////////////////////////////
324                     // The parameter to the Fjt function
325                     const SIMINT_DBLTYPE F_x = SIMINT_MUL(R2, alpha);
326 
327 
328                     const SIMINT_DBLTYPE Q_prefac = mask_load(nlane, Q.prefac + j);
329 
330 
331                     boys_F_split(PRIM_INT__s_s_s_s, F_x, 16);
332                     SIMINT_DBLTYPE prefac = SIMINT_SQRT(one_over_PQalpha_sum);
333                     prefac = SIMINT_MUL(SIMINT_MUL(P_prefac, Q_prefac), prefac);
334                     for(n = 0; n <= 16; n++)
335                         PRIM_INT__s_s_s_s[n] = SIMINT_MUL(PRIM_INT__s_s_s_s[n], prefac);
336 
337                     //////////////////////////////////////////////
338                     // Primitive integrals: Vertical recurrance
339                     //////////////////////////////////////////////
340 
341                     const SIMINT_DBLTYPE vrr_const_1_over_2p = one_over_2p;
342                     const SIMINT_DBLTYPE vrr_const_2_over_2p = SIMINT_MUL(const_2, one_over_2p);
343                     const SIMINT_DBLTYPE vrr_const_3_over_2p = SIMINT_MUL(const_3, one_over_2p);
344                     const SIMINT_DBLTYPE vrr_const_4_over_2p = SIMINT_MUL(const_4, one_over_2p);
345                     const SIMINT_DBLTYPE vrr_const_5_over_2p = SIMINT_MUL(const_5, one_over_2p);
346                     const SIMINT_DBLTYPE vrr_const_6_over_2p = SIMINT_MUL(const_6, one_over_2p);
347                     const SIMINT_DBLTYPE vrr_const_7_over_2p = SIMINT_MUL(const_7, one_over_2p);
348                     const SIMINT_DBLTYPE vrr_const_8_over_2p = SIMINT_MUL(const_8, one_over_2p);
349                     const SIMINT_DBLTYPE vrr_const_9_over_2p = SIMINT_MUL(const_9, one_over_2p);
350                     const SIMINT_DBLTYPE vrr_const_1_over_2q = one_over_2q;
351                     const SIMINT_DBLTYPE vrr_const_2_over_2q = SIMINT_MUL(const_2, one_over_2q);
352                     const SIMINT_DBLTYPE vrr_const_3_over_2q = SIMINT_MUL(const_3, one_over_2q);
353                     const SIMINT_DBLTYPE vrr_const_4_over_2q = SIMINT_MUL(const_4, one_over_2q);
354                     const SIMINT_DBLTYPE vrr_const_5_over_2q = SIMINT_MUL(const_5, one_over_2q);
355                     const SIMINT_DBLTYPE vrr_const_1_over_2pq = one_over_2pq;
356                     const SIMINT_DBLTYPE vrr_const_2_over_2pq = SIMINT_MUL(const_2, one_over_2pq);
357                     const SIMINT_DBLTYPE vrr_const_3_over_2pq = SIMINT_MUL(const_3, one_over_2pq);
358                     const SIMINT_DBLTYPE vrr_const_4_over_2pq = SIMINT_MUL(const_4, one_over_2pq);
359                     const SIMINT_DBLTYPE vrr_const_5_over_2pq = SIMINT_MUL(const_5, one_over_2pq);
360                     const SIMINT_DBLTYPE vrr_const_6_over_2pq = SIMINT_MUL(const_6, one_over_2pq);
361                     const SIMINT_DBLTYPE vrr_const_7_over_2pq = SIMINT_MUL(const_7, one_over_2pq);
362                     const SIMINT_DBLTYPE vrr_const_8_over_2pq = SIMINT_MUL(const_8, one_over_2pq);
363                     const SIMINT_DBLTYPE vrr_const_9_over_2pq = SIMINT_MUL(const_9, one_over_2pq);
364                     const SIMINT_DBLTYPE vrr_const_10_over_2pq = SIMINT_MUL(const_10, one_over_2pq);
365 
366 
367 
368                     // Forming PRIM_INT__p_s_s_s[16 * 3];
369                     for(n = 0; n < 16; ++n)  // loop over orders of auxiliary function
370                     {
371 
372                         PRIM_INT__p_s_s_s[n * 3 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
373                         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]);
374 
375                         PRIM_INT__p_s_s_s[n * 3 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
376                         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]);
377 
378                         PRIM_INT__p_s_s_s[n * 3 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
379                         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]);
380 
381                     }
382 
383 
384 
385                     // Forming PRIM_INT__d_s_s_s[15 * 6];
386                     for(n = 0; n < 15; ++n)  // loop over orders of auxiliary function
387                     {
388 
389                         PRIM_INT__d_s_s_s[n * 6 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__p_s_s_s[n * 3 + 0]);
390                         PRIM_INT__d_s_s_s[n * 6 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_s_s[n * 6 + 0]);
391                         PRIM_INT__d_s_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__d_s_s_s[n * 6 + 0]);
392 
393                         PRIM_INT__d_s_s_s[n * 6 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__p_s_s_s[n * 3 + 0]);
394                         PRIM_INT__d_s_s_s[n * 6 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_s_s[n * 6 + 1]);
395 
396                         PRIM_INT__d_s_s_s[n * 6 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 0]);
397                         PRIM_INT__d_s_s_s[n * 6 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_s_s[n * 6 + 2]);
398 
399                         PRIM_INT__d_s_s_s[n * 6 + 3] = SIMINT_MUL(P_PA[1], PRIM_INT__p_s_s_s[n * 3 + 1]);
400                         PRIM_INT__d_s_s_s[n * 6 + 3] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_s_s[n * 6 + 3]);
401                         PRIM_INT__d_s_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__d_s_s_s[n * 6 + 3]);
402 
403                         PRIM_INT__d_s_s_s[n * 6 + 4] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 1]);
404                         PRIM_INT__d_s_s_s[n * 6 + 4] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_s_s[n * 6 + 4]);
405 
406                         PRIM_INT__d_s_s_s[n * 6 + 5] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 2]);
407                         PRIM_INT__d_s_s_s[n * 6 + 5] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__d_s_s_s[n * 6 + 5]);
408                         PRIM_INT__d_s_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__d_s_s_s[n * 6 + 5]);
409 
410                     }
411 
412 
413 
414                     // Forming PRIM_INT__f_s_s_s[14 * 10];
415                     for(n = 0; n < 14; ++n)  // loop over orders of auxiliary function
416                     {
417 
418                         PRIM_INT__f_s_s_s[n * 10 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 0]);
419                         PRIM_INT__f_s_s_s[n * 10 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__f_s_s_s[n * 10 + 0]);
420                         PRIM_INT__f_s_s_s[n * 10 + 0] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_s_s[n * 3 + 0]), PRIM_INT__f_s_s_s[n * 10 + 0]);
421 
422                         PRIM_INT__f_s_s_s[n * 10 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 0]);
423                         PRIM_INT__f_s_s_s[n * 10 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__f_s_s_s[n * 10 + 1]);
424 
425                         PRIM_INT__f_s_s_s[n * 10 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 0]);
426                         PRIM_INT__f_s_s_s[n * 10 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__f_s_s_s[n * 10 + 2]);
427 
428                         PRIM_INT__f_s_s_s[n * 10 + 3] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 3]);
429                         PRIM_INT__f_s_s_s[n * 10 + 3] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__f_s_s_s[n * 10 + 3]);
430 
431                         PRIM_INT__f_s_s_s[n * 10 + 4] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 1]);
432                         PRIM_INT__f_s_s_s[n * 10 + 4] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__f_s_s_s[n * 10 + 4]);
433 
434                         PRIM_INT__f_s_s_s[n * 10 + 5] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 5]);
435                         PRIM_INT__f_s_s_s[n * 10 + 5] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__f_s_s_s[n * 10 + 5]);
436 
437                         PRIM_INT__f_s_s_s[n * 10 + 6] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 3]);
438                         PRIM_INT__f_s_s_s[n * 10 + 6] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__f_s_s_s[n * 10 + 6]);
439                         PRIM_INT__f_s_s_s[n * 10 + 6] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_s_s[n * 3 + 1]), PRIM_INT__f_s_s_s[n * 10 + 6]);
440 
441                         PRIM_INT__f_s_s_s[n * 10 + 7] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 3]);
442                         PRIM_INT__f_s_s_s[n * 10 + 7] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__f_s_s_s[n * 10 + 7]);
443 
444                         PRIM_INT__f_s_s_s[n * 10 + 8] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 5]);
445                         PRIM_INT__f_s_s_s[n * 10 + 8] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__f_s_s_s[n * 10 + 8]);
446 
447                         PRIM_INT__f_s_s_s[n * 10 + 9] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 5]);
448                         PRIM_INT__f_s_s_s[n * 10 + 9] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__f_s_s_s[n * 10 + 9]);
449                         PRIM_INT__f_s_s_s[n * 10 + 9] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_s_s[n * 3 + 2]), PRIM_INT__f_s_s_s[n * 10 + 9]);
450 
451                     }
452 
453 
454                     VRR_I_g_s_s_s(
455                             PRIM_INT__g_s_s_s,
456                             PRIM_INT__f_s_s_s,
457                             PRIM_INT__d_s_s_s,
458                             P_PA,
459                             a_over_p,
460                             aop_PQ,
461                             one_over_2p,
462                             13);
463 
464 
465                     VRR_I_h_s_s_s(
466                             PRIM_INT__h_s_s_s,
467                             PRIM_INT__g_s_s_s,
468                             PRIM_INT__f_s_s_s,
469                             P_PA,
470                             a_over_p,
471                             aop_PQ,
472                             one_over_2p,
473                             12);
474 
475 
476                     ostei_general_vrr_K(5, 0, 1, 0, 6,
477                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
478                             PRIM_INT__h_s_s_s, NULL, NULL, PRIM_INT__g_s_s_s, NULL, PRIM_INT__h_s_p_s);
479 
480 
481                     VRR_K_g_s_p_s(
482                             PRIM_INT__g_s_p_s,
483                             PRIM_INT__g_s_s_s,
484                             PRIM_INT__f_s_s_s,
485                             Q_PA,
486                             aoq_PQ,
487                             one_over_2pq,
488                             6);
489 
490 
491                     ostei_general_vrr_K(5, 0, 2, 0, 5,
492                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
493                             PRIM_INT__h_s_p_s, PRIM_INT__h_s_s_s, NULL, PRIM_INT__g_s_p_s, NULL, PRIM_INT__h_s_d_s);
494 
495 
496                     VRR_K_f_s_p_s(
497                             PRIM_INT__f_s_p_s,
498                             PRIM_INT__f_s_s_s,
499                             PRIM_INT__d_s_s_s,
500                             Q_PA,
501                             aoq_PQ,
502                             one_over_2pq,
503                             6);
504 
505 
506                     ostei_general_vrr_K(4, 0, 2, 0, 5,
507                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
508                             PRIM_INT__g_s_p_s, PRIM_INT__g_s_s_s, NULL, PRIM_INT__f_s_p_s, NULL, PRIM_INT__g_s_d_s);
509 
510 
511                     ostei_general_vrr_K(5, 0, 3, 0, 4,
512                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
513                             PRIM_INT__h_s_d_s, PRIM_INT__h_s_p_s, NULL, PRIM_INT__g_s_d_s, NULL, PRIM_INT__h_s_f_s);
514 
515 
516 
517                     // Forming PRIM_INT__d_s_p_s[6 * 18];
518                     for(n = 0; n < 6; ++n)  // loop over orders of auxiliary function
519                     {
520 
521                         PRIM_INT__d_s_p_s[n * 18 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 0]);
522                         PRIM_INT__d_s_p_s[n * 18 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__d_s_p_s[n * 18 + 0]);
523                         PRIM_INT__d_s_p_s[n * 18 + 0] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_p_s[n * 18 + 0]);
524 
525                         PRIM_INT__d_s_p_s[n * 18 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 0]);
526                         PRIM_INT__d_s_p_s[n * 18 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__d_s_p_s[n * 18 + 1]);
527 
528                         PRIM_INT__d_s_p_s[n * 18 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 0]);
529                         PRIM_INT__d_s_p_s[n * 18 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__d_s_p_s[n * 18 + 2]);
530 
531                         PRIM_INT__d_s_p_s[n * 18 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 1]);
532                         PRIM_INT__d_s_p_s[n * 18 + 3] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__d_s_p_s[n * 18 + 3]);
533                         PRIM_INT__d_s_p_s[n * 18 + 3] = 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 + 3]);
534 
535                         PRIM_INT__d_s_p_s[n * 18 + 4] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 1]);
536                         PRIM_INT__d_s_p_s[n * 18 + 4] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__d_s_p_s[n * 18 + 4]);
537                         PRIM_INT__d_s_p_s[n * 18 + 4] = 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 + 4]);
538 
539                         PRIM_INT__d_s_p_s[n * 18 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 1]);
540                         PRIM_INT__d_s_p_s[n * 18 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__d_s_p_s[n * 18 + 5]);
541 
542                         PRIM_INT__d_s_p_s[n * 18 + 6] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 2]);
543                         PRIM_INT__d_s_p_s[n * 18 + 6] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 2], PRIM_INT__d_s_p_s[n * 18 + 6]);
544                         PRIM_INT__d_s_p_s[n * 18 + 6] = 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 + 6]);
545 
546                         PRIM_INT__d_s_p_s[n * 18 + 7] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 2]);
547                         PRIM_INT__d_s_p_s[n * 18 + 7] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 2], PRIM_INT__d_s_p_s[n * 18 + 7]);
548 
549                         PRIM_INT__d_s_p_s[n * 18 + 8] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 2]);
550                         PRIM_INT__d_s_p_s[n * 18 + 8] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 2], PRIM_INT__d_s_p_s[n * 18 + 8]);
551                         PRIM_INT__d_s_p_s[n * 18 + 8] = 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 + 8]);
552 
553                         PRIM_INT__d_s_p_s[n * 18 + 9] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 3]);
554                         PRIM_INT__d_s_p_s[n * 18 + 9] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__d_s_p_s[n * 18 + 9]);
555 
556                         PRIM_INT__d_s_p_s[n * 18 + 10] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 3]);
557                         PRIM_INT__d_s_p_s[n * 18 + 10] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__d_s_p_s[n * 18 + 10]);
558                         PRIM_INT__d_s_p_s[n * 18 + 10] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_p_s[n * 18 + 10]);
559 
560                         PRIM_INT__d_s_p_s[n * 18 + 11] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 3]);
561                         PRIM_INT__d_s_p_s[n * 18 + 11] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__d_s_p_s[n * 18 + 11]);
562 
563                         PRIM_INT__d_s_p_s[n * 18 + 12] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 4]);
564                         PRIM_INT__d_s_p_s[n * 18 + 12] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 4], PRIM_INT__d_s_p_s[n * 18 + 12]);
565 
566                         PRIM_INT__d_s_p_s[n * 18 + 13] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 4]);
567                         PRIM_INT__d_s_p_s[n * 18 + 13] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 4], PRIM_INT__d_s_p_s[n * 18 + 13]);
568                         PRIM_INT__d_s_p_s[n * 18 + 13] = 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 + 13]);
569 
570                         PRIM_INT__d_s_p_s[n * 18 + 14] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 4]);
571                         PRIM_INT__d_s_p_s[n * 18 + 14] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 4], PRIM_INT__d_s_p_s[n * 18 + 14]);
572                         PRIM_INT__d_s_p_s[n * 18 + 14] = 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 + 14]);
573 
574                         PRIM_INT__d_s_p_s[n * 18 + 15] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 5]);
575                         PRIM_INT__d_s_p_s[n * 18 + 15] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__d_s_p_s[n * 18 + 15]);
576 
577                         PRIM_INT__d_s_p_s[n * 18 + 16] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 5]);
578                         PRIM_INT__d_s_p_s[n * 18 + 16] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__d_s_p_s[n * 18 + 16]);
579 
580                         PRIM_INT__d_s_p_s[n * 18 + 17] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 5]);
581                         PRIM_INT__d_s_p_s[n * 18 + 17] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__d_s_p_s[n * 18 + 17]);
582                         PRIM_INT__d_s_p_s[n * 18 + 17] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__d_s_p_s[n * 18 + 17]);
583 
584                     }
585 
586 
587                     VRR_K_f_s_d_s(
588                             PRIM_INT__f_s_d_s,
589                             PRIM_INT__f_s_p_s,
590                             PRIM_INT__f_s_s_s,
591                             PRIM_INT__d_s_p_s,
592                             Q_PA,
593                             a_over_q,
594                             aoq_PQ,
595                             one_over_2pq,
596                             one_over_2q,
597                             5);
598 
599 
600                     ostei_general_vrr_K(4, 0, 3, 0, 4,
601                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
602                             PRIM_INT__g_s_d_s, PRIM_INT__g_s_p_s, NULL, PRIM_INT__f_s_d_s, NULL, PRIM_INT__g_s_f_s);
603 
604 
605                     ostei_general_vrr_K(5, 0, 4, 0, 3,
606                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
607                             PRIM_INT__h_s_f_s, PRIM_INT__h_s_d_s, NULL, PRIM_INT__g_s_f_s, NULL, PRIM_INT__h_s_g_s);
608 
609 
610 
611                     // Forming PRIM_INT__p_s_p_s[6 * 9];
612                     for(n = 0; n < 6; ++n)  // loop over orders of auxiliary function
613                     {
614 
615                         PRIM_INT__p_s_p_s[n * 9 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_s_s[n * 3 + 0]);
616                         PRIM_INT__p_s_p_s[n * 9 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_p_s[n * 9 + 0]);
617                         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]);
618 
619                         PRIM_INT__p_s_p_s[n * 9 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_s_s[n * 3 + 0]);
620                         PRIM_INT__p_s_p_s[n * 9 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_p_s[n * 9 + 1]);
621 
622                         PRIM_INT__p_s_p_s[n * 9 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_s_s[n * 3 + 0]);
623                         PRIM_INT__p_s_p_s[n * 9 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_p_s[n * 9 + 2]);
624 
625                         PRIM_INT__p_s_p_s[n * 9 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_s_s[n * 3 + 1]);
626                         PRIM_INT__p_s_p_s[n * 9 + 3] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_p_s[n * 9 + 3]);
627 
628                         PRIM_INT__p_s_p_s[n * 9 + 4] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_s_s[n * 3 + 1]);
629                         PRIM_INT__p_s_p_s[n * 9 + 4] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_p_s[n * 9 + 4]);
630                         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]);
631 
632                         PRIM_INT__p_s_p_s[n * 9 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_s_s[n * 3 + 1]);
633                         PRIM_INT__p_s_p_s[n * 9 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_p_s[n * 9 + 5]);
634 
635                         PRIM_INT__p_s_p_s[n * 9 + 6] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_s_s[n * 3 + 2]);
636                         PRIM_INT__p_s_p_s[n * 9 + 6] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_p_s[n * 9 + 6]);
637 
638                         PRIM_INT__p_s_p_s[n * 9 + 7] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_s_s[n * 3 + 2]);
639                         PRIM_INT__p_s_p_s[n * 9 + 7] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_p_s[n * 9 + 7]);
640 
641                         PRIM_INT__p_s_p_s[n * 9 + 8] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_s_s[n * 3 + 2]);
642                         PRIM_INT__p_s_p_s[n * 9 + 8] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_p_s[n * 9 + 8]);
643                         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]);
644 
645                     }
646 
647 
648                     VRR_K_d_s_d_s(
649                             PRIM_INT__d_s_d_s,
650                             PRIM_INT__d_s_p_s,
651                             PRIM_INT__d_s_s_s,
652                             PRIM_INT__p_s_p_s,
653                             Q_PA,
654                             a_over_q,
655                             aoq_PQ,
656                             one_over_2pq,
657                             one_over_2q,
658                             5);
659 
660 
661                     ostei_general_vrr_K(3, 0, 3, 0, 4,
662                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
663                             PRIM_INT__f_s_d_s, PRIM_INT__f_s_p_s, NULL, PRIM_INT__d_s_d_s, NULL, PRIM_INT__f_s_f_s);
664 
665 
666                     ostei_general_vrr_K(4, 0, 4, 0, 3,
667                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
668                             PRIM_INT__g_s_f_s, PRIM_INT__g_s_d_s, NULL, PRIM_INT__f_s_f_s, NULL, PRIM_INT__g_s_g_s);
669 
670 
671                     ostei_general_vrr_K(5, 0, 5, 0, 2,
672                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
673                             PRIM_INT__h_s_g_s, PRIM_INT__h_s_f_s, NULL, PRIM_INT__g_s_g_s, NULL, PRIM_INT__h_s_h_s);
674 
675 
676 
677                     // Forming PRIM_INT__s_s_p_s[6 * 3];
678                     for(n = 0; n < 6; ++n)  // loop over orders of auxiliary function
679                     {
680 
681                         PRIM_INT__s_s_p_s[n * 3 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
682                         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]);
683 
684                         PRIM_INT__s_s_p_s[n * 3 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
685                         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]);
686 
687                         PRIM_INT__s_s_p_s[n * 3 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
688                         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]);
689 
690                     }
691 
692 
693 
694                     // Forming PRIM_INT__p_s_d_s[5 * 18];
695                     for(n = 0; n < 5; ++n)  // loop over orders of auxiliary function
696                     {
697 
698                         PRIM_INT__p_s_d_s[n * 18 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_p_s[n * 9 + 0]);
699                         PRIM_INT__p_s_d_s[n * 18 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__p_s_p_s[(n+1) * 9 + 0], PRIM_INT__p_s_d_s[n * 18 + 0]);
700                         PRIM_INT__p_s_d_s[n * 18 + 0] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_s_s[n * 3 + 0]), PRIM_INT__p_s_d_s[n * 18 + 0]);
701                         PRIM_INT__p_s_d_s[n * 18 + 0] = 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 + 0]);
702 
703                         PRIM_INT__p_s_d_s[n * 18 + 3] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_p_s[n * 9 + 1]);
704                         PRIM_INT__p_s_d_s[n * 18 + 3] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__p_s_p_s[(n+1) * 9 + 1], PRIM_INT__p_s_d_s[n * 18 + 3]);
705                         PRIM_INT__p_s_d_s[n * 18 + 3] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_s_s[n * 3 + 0]), PRIM_INT__p_s_d_s[n * 18 + 3]);
706 
707                         PRIM_INT__p_s_d_s[n * 18 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_p_s[n * 9 + 2]);
708                         PRIM_INT__p_s_d_s[n * 18 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__p_s_p_s[(n+1) * 9 + 2], PRIM_INT__p_s_d_s[n * 18 + 5]);
709                         PRIM_INT__p_s_d_s[n * 18 + 5] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_s_s[n * 3 + 0]), PRIM_INT__p_s_d_s[n * 18 + 5]);
710 
711                         PRIM_INT__p_s_d_s[n * 18 + 6] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_p_s[n * 9 + 3]);
712                         PRIM_INT__p_s_d_s[n * 18 + 6] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__p_s_p_s[(n+1) * 9 + 3], PRIM_INT__p_s_d_s[n * 18 + 6]);
713                         PRIM_INT__p_s_d_s[n * 18 + 6] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_s_s[n * 3 + 1]), PRIM_INT__p_s_d_s[n * 18 + 6]);
714 
715                         PRIM_INT__p_s_d_s[n * 18 + 9] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_p_s[n * 9 + 4]);
716                         PRIM_INT__p_s_d_s[n * 18 + 9] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__p_s_p_s[(n+1) * 9 + 4], PRIM_INT__p_s_d_s[n * 18 + 9]);
717                         PRIM_INT__p_s_d_s[n * 18 + 9] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_s_s[n * 3 + 1]), PRIM_INT__p_s_d_s[n * 18 + 9]);
718                         PRIM_INT__p_s_d_s[n * 18 + 9] = 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 + 9]);
719 
720                         PRIM_INT__p_s_d_s[n * 18 + 11] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_p_s[n * 9 + 5]);
721                         PRIM_INT__p_s_d_s[n * 18 + 11] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__p_s_p_s[(n+1) * 9 + 5], PRIM_INT__p_s_d_s[n * 18 + 11]);
722                         PRIM_INT__p_s_d_s[n * 18 + 11] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_s_s[n * 3 + 1]), PRIM_INT__p_s_d_s[n * 18 + 11]);
723 
724                         PRIM_INT__p_s_d_s[n * 18 + 12] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_p_s[n * 9 + 6]);
725                         PRIM_INT__p_s_d_s[n * 18 + 12] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__p_s_p_s[(n+1) * 9 + 6], PRIM_INT__p_s_d_s[n * 18 + 12]);
726                         PRIM_INT__p_s_d_s[n * 18 + 12] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_s_s[n * 3 + 2]), PRIM_INT__p_s_d_s[n * 18 + 12]);
727 
728                         PRIM_INT__p_s_d_s[n * 18 + 15] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_p_s[n * 9 + 7]);
729                         PRIM_INT__p_s_d_s[n * 18 + 15] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__p_s_p_s[(n+1) * 9 + 7], PRIM_INT__p_s_d_s[n * 18 + 15]);
730                         PRIM_INT__p_s_d_s[n * 18 + 15] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_s_s[n * 3 + 2]), PRIM_INT__p_s_d_s[n * 18 + 15]);
731 
732                         PRIM_INT__p_s_d_s[n * 18 + 17] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_p_s[n * 9 + 8]);
733                         PRIM_INT__p_s_d_s[n * 18 + 17] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__p_s_p_s[(n+1) * 9 + 8], PRIM_INT__p_s_d_s[n * 18 + 17]);
734                         PRIM_INT__p_s_d_s[n * 18 + 17] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_s_s[n * 3 + 2]), PRIM_INT__p_s_d_s[n * 18 + 17]);
735                         PRIM_INT__p_s_d_s[n * 18 + 17] = 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 + 17]);
736 
737                     }
738 
739 
740                     VRR_K_d_s_f_s(
741                             PRIM_INT__d_s_f_s,
742                             PRIM_INT__d_s_d_s,
743                             PRIM_INT__d_s_p_s,
744                             PRIM_INT__p_s_d_s,
745                             Q_PA,
746                             a_over_q,
747                             aoq_PQ,
748                             one_over_2pq,
749                             one_over_2q,
750                             4);
751 
752 
753                     ostei_general_vrr_K(3, 0, 4, 0, 3,
754                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
755                             PRIM_INT__f_s_f_s, PRIM_INT__f_s_d_s, NULL, PRIM_INT__d_s_f_s, NULL, PRIM_INT__f_s_g_s);
756 
757 
758                     ostei_general_vrr_K(4, 0, 5, 0, 2,
759                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
760                             PRIM_INT__g_s_g_s, PRIM_INT__g_s_f_s, NULL, PRIM_INT__f_s_g_s, NULL, PRIM_INT__g_s_h_s);
761 
762 
763                     ostei_general_vrr_K(5, 0, 6, 0, 1,
764                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
765                             PRIM_INT__h_s_h_s, PRIM_INT__h_s_g_s, NULL, PRIM_INT__g_s_h_s, NULL, PRIM_INT__h_s_i_s);
766 
767 
768                     ostei_general_vrr1_I(6, 11,
769                             one_over_2p, a_over_p, aop_PQ, P_PA,
770                             PRIM_INT__h_s_s_s, PRIM_INT__g_s_s_s, PRIM_INT__i_s_s_s);
771 
772 
773                     ostei_general_vrr_K(6, 0, 1, 0, 6,
774                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
775                             PRIM_INT__i_s_s_s, NULL, NULL, PRIM_INT__h_s_s_s, NULL, PRIM_INT__i_s_p_s);
776 
777 
778                     ostei_general_vrr_K(6, 0, 2, 0, 5,
779                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
780                             PRIM_INT__i_s_p_s, PRIM_INT__i_s_s_s, NULL, PRIM_INT__h_s_p_s, NULL, PRIM_INT__i_s_d_s);
781 
782 
783                     ostei_general_vrr_K(6, 0, 3, 0, 4,
784                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
785                             PRIM_INT__i_s_d_s, PRIM_INT__i_s_p_s, NULL, PRIM_INT__h_s_d_s, NULL, PRIM_INT__i_s_f_s);
786 
787 
788                     ostei_general_vrr_K(6, 0, 4, 0, 3,
789                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
790                             PRIM_INT__i_s_f_s, PRIM_INT__i_s_d_s, NULL, PRIM_INT__h_s_f_s, NULL, PRIM_INT__i_s_g_s);
791 
792 
793                     ostei_general_vrr_K(6, 0, 5, 0, 2,
794                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
795                             PRIM_INT__i_s_g_s, PRIM_INT__i_s_f_s, NULL, PRIM_INT__h_s_g_s, NULL, PRIM_INT__i_s_h_s);
796 
797 
798                     ostei_general_vrr_K(6, 0, 6, 0, 1,
799                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
800                             PRIM_INT__i_s_h_s, PRIM_INT__i_s_g_s, NULL, PRIM_INT__h_s_h_s, NULL, PRIM_INT__i_s_i_s);
801 
802 
803                     ostei_general_vrr1_I(7, 10,
804                             one_over_2p, a_over_p, aop_PQ, P_PA,
805                             PRIM_INT__i_s_s_s, PRIM_INT__h_s_s_s, PRIM_INT__k_s_s_s);
806 
807 
808                     ostei_general_vrr_K(7, 0, 1, 0, 6,
809                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
810                             PRIM_INT__k_s_s_s, NULL, NULL, PRIM_INT__i_s_s_s, NULL, PRIM_INT__k_s_p_s);
811 
812 
813                     ostei_general_vrr_K(7, 0, 2, 0, 5,
814                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
815                             PRIM_INT__k_s_p_s, PRIM_INT__k_s_s_s, NULL, PRIM_INT__i_s_p_s, NULL, PRIM_INT__k_s_d_s);
816 
817 
818                     ostei_general_vrr_K(7, 0, 3, 0, 4,
819                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
820                             PRIM_INT__k_s_d_s, PRIM_INT__k_s_p_s, NULL, PRIM_INT__i_s_d_s, NULL, PRIM_INT__k_s_f_s);
821 
822 
823                     ostei_general_vrr_K(7, 0, 4, 0, 3,
824                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
825                             PRIM_INT__k_s_f_s, PRIM_INT__k_s_d_s, NULL, PRIM_INT__i_s_f_s, NULL, PRIM_INT__k_s_g_s);
826 
827 
828                     ostei_general_vrr_K(7, 0, 5, 0, 2,
829                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
830                             PRIM_INT__k_s_g_s, PRIM_INT__k_s_f_s, NULL, PRIM_INT__i_s_g_s, NULL, PRIM_INT__k_s_h_s);
831 
832 
833                     ostei_general_vrr_K(7, 0, 6, 0, 1,
834                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
835                             PRIM_INT__k_s_h_s, PRIM_INT__k_s_g_s, NULL, PRIM_INT__i_s_h_s, NULL, PRIM_INT__k_s_i_s);
836 
837 
838                     ostei_general_vrr1_I(8, 9,
839                             one_over_2p, a_over_p, aop_PQ, P_PA,
840                             PRIM_INT__k_s_s_s, PRIM_INT__i_s_s_s, PRIM_INT__l_s_s_s);
841 
842 
843                     ostei_general_vrr_K(8, 0, 1, 0, 6,
844                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
845                             PRIM_INT__l_s_s_s, NULL, NULL, PRIM_INT__k_s_s_s, NULL, PRIM_INT__l_s_p_s);
846 
847 
848                     ostei_general_vrr_K(8, 0, 2, 0, 5,
849                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
850                             PRIM_INT__l_s_p_s, PRIM_INT__l_s_s_s, NULL, PRIM_INT__k_s_p_s, NULL, PRIM_INT__l_s_d_s);
851 
852 
853                     ostei_general_vrr_K(8, 0, 3, 0, 4,
854                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
855                             PRIM_INT__l_s_d_s, PRIM_INT__l_s_p_s, NULL, PRIM_INT__k_s_d_s, NULL, PRIM_INT__l_s_f_s);
856 
857 
858                     ostei_general_vrr_K(8, 0, 4, 0, 3,
859                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
860                             PRIM_INT__l_s_f_s, PRIM_INT__l_s_d_s, NULL, PRIM_INT__k_s_f_s, NULL, PRIM_INT__l_s_g_s);
861 
862 
863                     ostei_general_vrr_K(8, 0, 5, 0, 2,
864                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
865                             PRIM_INT__l_s_g_s, PRIM_INT__l_s_f_s, NULL, PRIM_INT__k_s_g_s, NULL, PRIM_INT__l_s_h_s);
866 
867 
868                     ostei_general_vrr_K(8, 0, 6, 0, 1,
869                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
870                             PRIM_INT__l_s_h_s, PRIM_INT__l_s_g_s, NULL, PRIM_INT__k_s_h_s, NULL, PRIM_INT__l_s_i_s);
871 
872 
873                     ostei_general_vrr1_I(9, 8,
874                             one_over_2p, a_over_p, aop_PQ, P_PA,
875                             PRIM_INT__l_s_s_s, PRIM_INT__k_s_s_s, PRIM_INT__m_s_s_s);
876 
877 
878                     ostei_general_vrr_K(9, 0, 1, 0, 6,
879                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
880                             PRIM_INT__m_s_s_s, NULL, NULL, PRIM_INT__l_s_s_s, NULL, PRIM_INT__m_s_p_s);
881 
882 
883                     ostei_general_vrr_K(9, 0, 2, 0, 5,
884                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
885                             PRIM_INT__m_s_p_s, PRIM_INT__m_s_s_s, NULL, PRIM_INT__l_s_p_s, NULL, PRIM_INT__m_s_d_s);
886 
887 
888                     ostei_general_vrr_K(9, 0, 3, 0, 4,
889                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
890                             PRIM_INT__m_s_d_s, PRIM_INT__m_s_p_s, NULL, PRIM_INT__l_s_d_s, NULL, PRIM_INT__m_s_f_s);
891 
892 
893                     ostei_general_vrr_K(9, 0, 4, 0, 3,
894                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
895                             PRIM_INT__m_s_f_s, PRIM_INT__m_s_d_s, NULL, PRIM_INT__l_s_f_s, NULL, PRIM_INT__m_s_g_s);
896 
897 
898                     ostei_general_vrr_K(9, 0, 5, 0, 2,
899                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
900                             PRIM_INT__m_s_g_s, PRIM_INT__m_s_f_s, NULL, PRIM_INT__l_s_g_s, NULL, PRIM_INT__m_s_h_s);
901 
902 
903                     ostei_general_vrr_K(9, 0, 6, 0, 1,
904                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
905                             PRIM_INT__m_s_h_s, PRIM_INT__m_s_g_s, NULL, PRIM_INT__l_s_h_s, NULL, PRIM_INT__m_s_i_s);
906 
907 
908                     ostei_general_vrr1_I(10, 7,
909                             one_over_2p, a_over_p, aop_PQ, P_PA,
910                             PRIM_INT__m_s_s_s, PRIM_INT__l_s_s_s, PRIM_INT__n_s_s_s);
911 
912 
913                     ostei_general_vrr_K(10, 0, 1, 0, 6,
914                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
915                             PRIM_INT__n_s_s_s, NULL, NULL, PRIM_INT__m_s_s_s, NULL, PRIM_INT__n_s_p_s);
916 
917 
918                     ostei_general_vrr_K(10, 0, 2, 0, 5,
919                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
920                             PRIM_INT__n_s_p_s, PRIM_INT__n_s_s_s, NULL, PRIM_INT__m_s_p_s, NULL, PRIM_INT__n_s_d_s);
921 
922 
923                     ostei_general_vrr_K(10, 0, 3, 0, 4,
924                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
925                             PRIM_INT__n_s_d_s, PRIM_INT__n_s_p_s, NULL, PRIM_INT__m_s_d_s, NULL, PRIM_INT__n_s_f_s);
926 
927 
928                     ostei_general_vrr_K(10, 0, 4, 0, 3,
929                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
930                             PRIM_INT__n_s_f_s, PRIM_INT__n_s_d_s, NULL, PRIM_INT__m_s_f_s, NULL, PRIM_INT__n_s_g_s);
931 
932 
933                     ostei_general_vrr_K(10, 0, 5, 0, 2,
934                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
935                             PRIM_INT__n_s_g_s, PRIM_INT__n_s_f_s, NULL, PRIM_INT__m_s_g_s, NULL, PRIM_INT__n_s_h_s);
936 
937 
938                     ostei_general_vrr_K(10, 0, 6, 0, 1,
939                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
940                             PRIM_INT__n_s_h_s, PRIM_INT__n_s_g_s, NULL, PRIM_INT__m_s_h_s, NULL, PRIM_INT__n_s_i_s);
941 
942 
943 
944 
945                     ////////////////////////////////////
946                     // Accumulate contracted integrals
947                     ////////////////////////////////////
948                     if(lastoffset == 0)
949                     {
950                         contract_all(441, PRIM_INT__h_s_h_s, PRIM_PTR_INT__h_s_h_s);
951                         contract_all(588, PRIM_INT__h_s_i_s, PRIM_PTR_INT__h_s_i_s);
952                         contract_all(588, PRIM_INT__i_s_h_s, PRIM_PTR_INT__i_s_h_s);
953                         contract_all(784, PRIM_INT__i_s_i_s, PRIM_PTR_INT__i_s_i_s);
954                         contract_all(756, PRIM_INT__k_s_h_s, PRIM_PTR_INT__k_s_h_s);
955                         contract_all(1008, PRIM_INT__k_s_i_s, PRIM_PTR_INT__k_s_i_s);
956                         contract_all(945, PRIM_INT__l_s_h_s, PRIM_PTR_INT__l_s_h_s);
957                         contract_all(1260, PRIM_INT__l_s_i_s, PRIM_PTR_INT__l_s_i_s);
958                         contract_all(1155, PRIM_INT__m_s_h_s, PRIM_PTR_INT__m_s_h_s);
959                         contract_all(1540, PRIM_INT__m_s_i_s, PRIM_PTR_INT__m_s_i_s);
960                         contract_all(1386, PRIM_INT__n_s_h_s, PRIM_PTR_INT__n_s_h_s);
961                         contract_all(1848, PRIM_INT__n_s_i_s, PRIM_PTR_INT__n_s_i_s);
962                     }
963                     else
964                     {
965                         contract(441, shelloffsets, PRIM_INT__h_s_h_s, PRIM_PTR_INT__h_s_h_s);
966                         contract(588, shelloffsets, PRIM_INT__h_s_i_s, PRIM_PTR_INT__h_s_i_s);
967                         contract(588, shelloffsets, PRIM_INT__i_s_h_s, PRIM_PTR_INT__i_s_h_s);
968                         contract(784, shelloffsets, PRIM_INT__i_s_i_s, PRIM_PTR_INT__i_s_i_s);
969                         contract(756, shelloffsets, PRIM_INT__k_s_h_s, PRIM_PTR_INT__k_s_h_s);
970                         contract(1008, shelloffsets, PRIM_INT__k_s_i_s, PRIM_PTR_INT__k_s_i_s);
971                         contract(945, shelloffsets, PRIM_INT__l_s_h_s, PRIM_PTR_INT__l_s_h_s);
972                         contract(1260, shelloffsets, PRIM_INT__l_s_i_s, PRIM_PTR_INT__l_s_i_s);
973                         contract(1155, shelloffsets, PRIM_INT__m_s_h_s, PRIM_PTR_INT__m_s_h_s);
974                         contract(1540, shelloffsets, PRIM_INT__m_s_i_s, PRIM_PTR_INT__m_s_i_s);
975                         contract(1386, shelloffsets, PRIM_INT__n_s_h_s, PRIM_PTR_INT__n_s_h_s);
976                         contract(1848, shelloffsets, PRIM_INT__n_s_i_s, PRIM_PTR_INT__n_s_i_s);
977                         PRIM_PTR_INT__h_s_h_s += lastoffset*441;
978                         PRIM_PTR_INT__h_s_i_s += lastoffset*588;
979                         PRIM_PTR_INT__i_s_h_s += lastoffset*588;
980                         PRIM_PTR_INT__i_s_i_s += lastoffset*784;
981                         PRIM_PTR_INT__k_s_h_s += lastoffset*756;
982                         PRIM_PTR_INT__k_s_i_s += lastoffset*1008;
983                         PRIM_PTR_INT__l_s_h_s += lastoffset*945;
984                         PRIM_PTR_INT__l_s_i_s += lastoffset*1260;
985                         PRIM_PTR_INT__m_s_h_s += lastoffset*1155;
986                         PRIM_PTR_INT__m_s_i_s += lastoffset*1540;
987                         PRIM_PTR_INT__n_s_h_s += lastoffset*1386;
988                         PRIM_PTR_INT__n_s_i_s += lastoffset*1848;
989                     }
990 
991                 }  // close loop over j
992             }  // close loop over i
993 
994             //Advance to the next batch
995             jstart = SIMINT_SIMD_ROUND(jend);
996 
997             //////////////////////////////////////////////
998             // Contracted integrals: Horizontal recurrance
999             //////////////////////////////////////////////
1000 
1001 
1002             const double hAB[3] = { P.AB_x[ab], P.AB_y[ab], P.AB_z[ab] };
1003 
1004 
1005             for(abcd = 0; abcd < nshellbatch; ++abcd, ++real_abcd)
1006             {
1007                 const double hCD[3] = { Q.AB_x[cd+abcd], Q.AB_y[cd+abcd], Q.AB_z[cd+abcd] };
1008 
1009                 // set up HRR pointers
1010                 double const * restrict HRR_INT__h_s_h_s = INT__h_s_h_s + abcd * 441;
1011                 double const * restrict HRR_INT__h_s_i_s = INT__h_s_i_s + abcd * 588;
1012                 double const * restrict HRR_INT__i_s_h_s = INT__i_s_h_s + abcd * 588;
1013                 double const * restrict HRR_INT__i_s_i_s = INT__i_s_i_s + abcd * 784;
1014                 double const * restrict HRR_INT__k_s_h_s = INT__k_s_h_s + abcd * 756;
1015                 double const * restrict HRR_INT__k_s_i_s = INT__k_s_i_s + abcd * 1008;
1016                 double const * restrict HRR_INT__l_s_h_s = INT__l_s_h_s + abcd * 945;
1017                 double const * restrict HRR_INT__l_s_i_s = INT__l_s_i_s + abcd * 1260;
1018                 double const * restrict HRR_INT__m_s_h_s = INT__m_s_h_s + abcd * 1155;
1019                 double const * restrict HRR_INT__m_s_i_s = INT__m_s_i_s + abcd * 1540;
1020                 double const * restrict HRR_INT__n_s_h_s = INT__n_s_h_s + abcd * 1386;
1021                 double const * restrict HRR_INT__n_s_i_s = INT__n_s_i_s + abcd * 1848;
1022                 double * restrict HRR_INT__h_h_h_p = INT__h_h_h_p + real_abcd * 27783;
1023 
1024                 // form INT__h_p_h_s
1025                 ostei_general_hrr_J(5, 1, 5, 0, hAB, HRR_INT__i_s_h_s, HRR_INT__h_s_h_s, HRR_INT__h_p_h_s);
1026 
1027                 // form INT__h_p_i_s
1028                 ostei_general_hrr_J(5, 1, 6, 0, hAB, HRR_INT__i_s_i_s, HRR_INT__h_s_i_s, HRR_INT__h_p_i_s);
1029 
1030                 // form INT__i_p_h_s
1031                 ostei_general_hrr_J(6, 1, 5, 0, hAB, HRR_INT__k_s_h_s, HRR_INT__i_s_h_s, HRR_INT__i_p_h_s);
1032 
1033                 // form INT__i_p_i_s
1034                 ostei_general_hrr_J(6, 1, 6, 0, hAB, HRR_INT__k_s_i_s, HRR_INT__i_s_i_s, HRR_INT__i_p_i_s);
1035 
1036                 // form INT__k_p_h_s
1037                 ostei_general_hrr_J(7, 1, 5, 0, hAB, HRR_INT__l_s_h_s, HRR_INT__k_s_h_s, HRR_INT__k_p_h_s);
1038 
1039                 // form INT__k_p_i_s
1040                 ostei_general_hrr_J(7, 1, 6, 0, hAB, HRR_INT__l_s_i_s, HRR_INT__k_s_i_s, HRR_INT__k_p_i_s);
1041 
1042                 // form INT__l_p_h_s
1043                 ostei_general_hrr_J(8, 1, 5, 0, hAB, HRR_INT__m_s_h_s, HRR_INT__l_s_h_s, HRR_INT__l_p_h_s);
1044 
1045                 // form INT__l_p_i_s
1046                 ostei_general_hrr_J(8, 1, 6, 0, hAB, HRR_INT__m_s_i_s, HRR_INT__l_s_i_s, HRR_INT__l_p_i_s);
1047 
1048                 // form INT__m_p_h_s
1049                 ostei_general_hrr_J(9, 1, 5, 0, hAB, HRR_INT__n_s_h_s, HRR_INT__m_s_h_s, HRR_INT__m_p_h_s);
1050 
1051                 // form INT__m_p_i_s
1052                 ostei_general_hrr_J(9, 1, 6, 0, hAB, HRR_INT__n_s_i_s, HRR_INT__m_s_i_s, HRR_INT__m_p_i_s);
1053 
1054                 // form INT__h_d_h_s
1055                 ostei_general_hrr_J(5, 2, 5, 0, hAB, HRR_INT__i_p_h_s, HRR_INT__h_p_h_s, HRR_INT__h_d_h_s);
1056 
1057                 // form INT__h_d_i_s
1058                 ostei_general_hrr_J(5, 2, 6, 0, hAB, HRR_INT__i_p_i_s, HRR_INT__h_p_i_s, HRR_INT__h_d_i_s);
1059 
1060                 // form INT__i_d_h_s
1061                 ostei_general_hrr_J(6, 2, 5, 0, hAB, HRR_INT__k_p_h_s, HRR_INT__i_p_h_s, HRR_INT__i_d_h_s);
1062 
1063                 // form INT__i_d_i_s
1064                 ostei_general_hrr_J(6, 2, 6, 0, hAB, HRR_INT__k_p_i_s, HRR_INT__i_p_i_s, HRR_INT__i_d_i_s);
1065 
1066                 // form INT__k_d_h_s
1067                 ostei_general_hrr_J(7, 2, 5, 0, hAB, HRR_INT__l_p_h_s, HRR_INT__k_p_h_s, HRR_INT__k_d_h_s);
1068 
1069                 // form INT__k_d_i_s
1070                 ostei_general_hrr_J(7, 2, 6, 0, hAB, HRR_INT__l_p_i_s, HRR_INT__k_p_i_s, HRR_INT__k_d_i_s);
1071 
1072                 // form INT__l_d_h_s
1073                 ostei_general_hrr_J(8, 2, 5, 0, hAB, HRR_INT__m_p_h_s, HRR_INT__l_p_h_s, HRR_INT__l_d_h_s);
1074 
1075                 // form INT__l_d_i_s
1076                 ostei_general_hrr_J(8, 2, 6, 0, hAB, HRR_INT__m_p_i_s, HRR_INT__l_p_i_s, HRR_INT__l_d_i_s);
1077 
1078                 // form INT__h_f_h_s
1079                 ostei_general_hrr_J(5, 3, 5, 0, hAB, HRR_INT__i_d_h_s, HRR_INT__h_d_h_s, HRR_INT__h_f_h_s);
1080 
1081                 // form INT__h_f_i_s
1082                 ostei_general_hrr_J(5, 3, 6, 0, hAB, HRR_INT__i_d_i_s, HRR_INT__h_d_i_s, HRR_INT__h_f_i_s);
1083 
1084                 // form INT__i_f_h_s
1085                 ostei_general_hrr_J(6, 3, 5, 0, hAB, HRR_INT__k_d_h_s, HRR_INT__i_d_h_s, HRR_INT__i_f_h_s);
1086 
1087                 // form INT__i_f_i_s
1088                 ostei_general_hrr_J(6, 3, 6, 0, hAB, HRR_INT__k_d_i_s, HRR_INT__i_d_i_s, HRR_INT__i_f_i_s);
1089 
1090                 // form INT__k_f_h_s
1091                 ostei_general_hrr_J(7, 3, 5, 0, hAB, HRR_INT__l_d_h_s, HRR_INT__k_d_h_s, HRR_INT__k_f_h_s);
1092 
1093                 // form INT__k_f_i_s
1094                 ostei_general_hrr_J(7, 3, 6, 0, hAB, HRR_INT__l_d_i_s, HRR_INT__k_d_i_s, HRR_INT__k_f_i_s);
1095 
1096                 // form INT__h_g_h_s
1097                 ostei_general_hrr_J(5, 4, 5, 0, hAB, HRR_INT__i_f_h_s, HRR_INT__h_f_h_s, HRR_INT__h_g_h_s);
1098 
1099                 // form INT__h_g_i_s
1100                 ostei_general_hrr_J(5, 4, 6, 0, hAB, HRR_INT__i_f_i_s, HRR_INT__h_f_i_s, HRR_INT__h_g_i_s);
1101 
1102                 // form INT__i_g_h_s
1103                 ostei_general_hrr_J(6, 4, 5, 0, hAB, HRR_INT__k_f_h_s, HRR_INT__i_f_h_s, HRR_INT__i_g_h_s);
1104 
1105                 // form INT__i_g_i_s
1106                 ostei_general_hrr_J(6, 4, 6, 0, hAB, HRR_INT__k_f_i_s, HRR_INT__i_f_i_s, HRR_INT__i_g_i_s);
1107 
1108                 // form INT__h_h_h_s
1109                 ostei_general_hrr_J(5, 5, 5, 0, hAB, HRR_INT__i_g_h_s, HRR_INT__h_g_h_s, HRR_INT__h_h_h_s);
1110 
1111                 // form INT__h_h_i_s
1112                 ostei_general_hrr_J(5, 5, 6, 0, hAB, HRR_INT__i_g_i_s, HRR_INT__h_g_i_s, HRR_INT__h_h_i_s);
1113 
1114                 // form INT__h_h_h_p
1115                 ostei_general_hrr_L(5, 5, 5, 1, hCD, HRR_INT__h_h_i_s, HRR_INT__h_h_h_s, HRR_INT__h_h_h_p);
1116 
1117 
1118             }  // close HRR loop
1119 
1120 
1121         }   // close loop cdbatch
1122 
1123         istart = iend;
1124     }  // close loop over ab
1125 
1126     return P.nshell12_clip * Q.nshell12_clip;
1127 }
1128 
ostei_h_h_p_h(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__h_h_p_h)1129 int ostei_h_h_p_h(struct simint_multi_shellpair const P,
1130                   struct simint_multi_shellpair const Q,
1131                   double screen_tol,
1132                   double * const restrict work,
1133                   double * const restrict INT__h_h_p_h)
1134 {
1135     double Q_AB[3*Q.nshell12];
1136     struct simint_multi_shellpair Q_tmp = Q;
1137     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
1138     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
1139     Q_tmp.AB_x = Q_AB;
1140     Q_tmp.AB_y = Q_AB + Q.nshell12;
1141     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
1142 
1143     for(int i = 0; i < Q.nshell12; i++)
1144     {
1145         Q_tmp.AB_x[i] = -Q.AB_x[i];
1146         Q_tmp.AB_y[i] = -Q.AB_y[i];
1147         Q_tmp.AB_z[i] = -Q.AB_z[i];
1148     }
1149 
1150     int ret = ostei_h_h_h_p(P, Q_tmp, screen_tol, work, INT__h_h_p_h);
1151     double buffer[27783] SIMINT_ALIGN_ARRAY_DBL;
1152 
1153     for(int q = 0; q < ret; q++)
1154     {
1155         int idx = 0;
1156         for(int a = 0; a < 21; ++a)
1157         for(int b = 0; b < 21; ++b)
1158         for(int c = 0; c < 3; ++c)
1159         for(int d = 0; d < 21; ++d)
1160             buffer[idx++] = INT__h_h_p_h[q*27783+a*1323+b*63+d*3+c];
1161 
1162         memcpy(INT__h_h_p_h+q*27783, buffer, 27783*sizeof(double));
1163     }
1164 
1165     return ret;
1166 }
1167 
1168