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_g_f_g_f(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__g_f_g_f)8 int ostei_g_f_g_f(struct simint_multi_shellpair const P,
9                   struct simint_multi_shellpair const Q,
10                   double screen_tol,
11                   double * const restrict work,
12                   double * const restrict INT__g_f_g_f)
13 {
14 
15     SIMINT_ASSUME_ALIGN_DBL(work);
16     SIMINT_ASSUME_ALIGN_DBL(INT__g_f_g_f);
17     int ab, cd, abcd;
18     int istart, jstart;
19     int iprimcd, nprim_icd, icd;
20     const int check_screen = (screen_tol > 0.0);
21     int i, j;
22     int n;
23     int not_screened;
24     int real_abcd;
25     int iket;
26     int ibra;
27 
28     // partition workspace
29     double * const INT__g_s_g_s = work + (SIMINT_NSHELL_SIMD * 0);
30     double * const INT__g_s_h_s = work + (SIMINT_NSHELL_SIMD * 225);
31     double * const INT__g_s_i_s = work + (SIMINT_NSHELL_SIMD * 540);
32     double * const INT__g_s_k_s = work + (SIMINT_NSHELL_SIMD * 960);
33     double * const INT__h_s_g_s = work + (SIMINT_NSHELL_SIMD * 1500);
34     double * const INT__h_s_h_s = work + (SIMINT_NSHELL_SIMD * 1815);
35     double * const INT__h_s_i_s = work + (SIMINT_NSHELL_SIMD * 2256);
36     double * const INT__h_s_k_s = work + (SIMINT_NSHELL_SIMD * 2844);
37     double * const INT__i_s_g_s = work + (SIMINT_NSHELL_SIMD * 3600);
38     double * const INT__i_s_h_s = work + (SIMINT_NSHELL_SIMD * 4020);
39     double * const INT__i_s_i_s = work + (SIMINT_NSHELL_SIMD * 4608);
40     double * const INT__i_s_k_s = work + (SIMINT_NSHELL_SIMD * 5392);
41     double * const INT__k_s_g_s = work + (SIMINT_NSHELL_SIMD * 6400);
42     double * const INT__k_s_h_s = work + (SIMINT_NSHELL_SIMD * 6940);
43     double * const INT__k_s_i_s = work + (SIMINT_NSHELL_SIMD * 7696);
44     double * const INT__k_s_k_s = work + (SIMINT_NSHELL_SIMD * 8704);
45     SIMINT_DBLTYPE * const primwork = (SIMINT_DBLTYPE *)(work + SIMINT_NSHELL_SIMD*10000);
46     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_s = primwork + 0;
47     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_p_s = primwork + 15;
48     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_d_s = primwork + 36;
49     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_f_s = primwork + 72;
50     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_s_s = primwork + 122;
51     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_p_s = primwork + 164;
52     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_d_s = primwork + 227;
53     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_f_s = primwork + 335;
54     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_g_s = primwork + 485;
55     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_s_s = primwork + 665;
56     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_p_s = primwork + 743;
57     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_d_s = primwork + 869;
58     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_f_s = primwork + 1085;
59     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_g_s = primwork + 1385;
60     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_h_s = primwork + 1745;
61     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_s_s = primwork + 2123;
62     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_p_s = primwork + 2243;
63     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_d_s = primwork + 2453;
64     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_f_s = primwork + 2813;
65     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_g_s = primwork + 3313;
66     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_h_s = primwork + 3913;
67     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_i_s = primwork + 4543;
68     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_s_s = primwork + 5103;
69     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_p_s = primwork + 5268;
70     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_d_s = primwork + 5583;
71     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_f_s = primwork + 6123;
72     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_g_s = primwork + 6873;
73     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_h_s = primwork + 7773;
74     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_i_s = primwork + 8718;
75     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_k_s = primwork + 9558;
76     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_s_s = primwork + 10098;
77     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_p_s = primwork + 10308;
78     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_d_s = primwork + 10749;
79     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_f_s = primwork + 11505;
80     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_g_s = primwork + 12555;
81     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_h_s = primwork + 13815;
82     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_i_s = primwork + 15138;
83     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_k_s = primwork + 16314;
84     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_s_s = primwork + 17070;
85     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_p_s = primwork + 17322;
86     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_d_s = primwork + 17910;
87     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_f_s = primwork + 18918;
88     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_g_s = primwork + 20318;
89     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_h_s = primwork + 21998;
90     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_i_s = primwork + 23762;
91     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_k_s = primwork + 25330;
92     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_s_s = primwork + 26338;
93     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_p_s = primwork + 26626;
94     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_d_s = primwork + 27382;
95     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_f_s = primwork + 28678;
96     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_g_s = primwork + 30478;
97     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_h_s = primwork + 32638;
98     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_i_s = primwork + 34906;
99     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_k_s = primwork + 36922;
100     double * const hrrwork = (double *)(primwork + 38218);
101     double * const HRR_INT__g_p_g_s = hrrwork + 0;
102     double * const HRR_INT__g_p_h_s = hrrwork + 675;
103     double * const HRR_INT__g_p_i_s = hrrwork + 1620;
104     double * const HRR_INT__g_p_k_s = hrrwork + 2880;
105     double * const HRR_INT__g_d_g_s = hrrwork + 4500;
106     double * const HRR_INT__g_d_h_s = hrrwork + 5850;
107     double * const HRR_INT__g_d_i_s = hrrwork + 7740;
108     double * const HRR_INT__g_d_k_s = hrrwork + 10260;
109     double * const HRR_INT__g_f_g_s = hrrwork + 13500;
110     double * const HRR_INT__g_f_g_p = hrrwork + 15750;
111     double * const HRR_INT__g_f_g_d = hrrwork + 22500;
112     double * const HRR_INT__g_f_h_s = hrrwork + 36000;
113     double * const HRR_INT__g_f_h_p = hrrwork + 39150;
114     double * const HRR_INT__g_f_h_d = hrrwork + 48600;
115     double * const HRR_INT__g_f_i_s = hrrwork + 67500;
116     double * const HRR_INT__g_f_i_p = hrrwork + 71700;
117     double * const HRR_INT__g_f_k_s = hrrwork + 84300;
118     double * const HRR_INT__h_p_g_s = hrrwork + 89700;
119     double * const HRR_INT__h_p_h_s = hrrwork + 90645;
120     double * const HRR_INT__h_p_i_s = hrrwork + 91968;
121     double * const HRR_INT__h_p_k_s = hrrwork + 93732;
122     double * const HRR_INT__h_d_g_s = hrrwork + 96000;
123     double * const HRR_INT__h_d_h_s = hrrwork + 97890;
124     double * const HRR_INT__h_d_i_s = hrrwork + 100536;
125     double * const HRR_INT__h_d_k_s = hrrwork + 104064;
126     double * const HRR_INT__i_p_g_s = hrrwork + 108600;
127     double * const HRR_INT__i_p_h_s = hrrwork + 109860;
128     double * const HRR_INT__i_p_i_s = hrrwork + 111624;
129     double * const HRR_INT__i_p_k_s = hrrwork + 113976;
130 
131 
132     // Create constants
133     const SIMINT_DBLTYPE const_1 = SIMINT_DBLSET1(1);
134     const SIMINT_DBLTYPE const_2 = SIMINT_DBLSET1(2);
135     const SIMINT_DBLTYPE const_3 = SIMINT_DBLSET1(3);
136     const SIMINT_DBLTYPE const_4 = SIMINT_DBLSET1(4);
137     const SIMINT_DBLTYPE const_5 = SIMINT_DBLSET1(5);
138     const SIMINT_DBLTYPE const_6 = SIMINT_DBLSET1(6);
139     const SIMINT_DBLTYPE const_7 = SIMINT_DBLSET1(7);
140     const SIMINT_DBLTYPE one_half = SIMINT_DBLSET1(0.5);
141 
142 
143     ////////////////////////////////////////
144     // Loop over shells and primitives
145     ////////////////////////////////////////
146 
147     real_abcd = 0;
148     istart = 0;
149     for(ab = 0; ab < P.nshell12_clip; ++ab)
150     {
151         const int iend = istart + P.nprim12[ab];
152 
153         cd = 0;
154         jstart = 0;
155 
156         for(cd = 0; cd < Q.nshell12_clip; cd += SIMINT_NSHELL_SIMD)
157         {
158             const int nshellbatch = ((cd + SIMINT_NSHELL_SIMD) > Q.nshell12_clip) ? Q.nshell12_clip - cd : SIMINT_NSHELL_SIMD;
159             int jend = jstart;
160             for(i = 0; i < nshellbatch; i++)
161                 jend += Q.nprim12[cd+i];
162 
163             // Clear the beginning of the workspace (where we are accumulating integrals)
164             memset(work, 0, SIMINT_NSHELL_SIMD * 10000 * sizeof(double));
165             abcd = 0;
166 
167 
168             for(i = istart; i < iend; ++i)
169             {
170                 SIMINT_DBLTYPE bra_screen_max;  // only used if check_screen
171 
172                 if(check_screen)
173                 {
174                     // Skip this whole thing if always insignificant
175                     if((P.screen[i] * Q.screen_max) < screen_tol)
176                         continue;
177                     bra_screen_max = SIMINT_DBLSET1(P.screen[i]);
178                 }
179 
180                 icd = 0;
181                 iprimcd = 0;
182                 nprim_icd = Q.nprim12[cd];
183                 double * restrict PRIM_PTR_INT__g_s_g_s = INT__g_s_g_s + abcd * 225;
184                 double * restrict PRIM_PTR_INT__g_s_h_s = INT__g_s_h_s + abcd * 315;
185                 double * restrict PRIM_PTR_INT__g_s_i_s = INT__g_s_i_s + abcd * 420;
186                 double * restrict PRIM_PTR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
187                 double * restrict PRIM_PTR_INT__h_s_g_s = INT__h_s_g_s + abcd * 315;
188                 double * restrict PRIM_PTR_INT__h_s_h_s = INT__h_s_h_s + abcd * 441;
189                 double * restrict PRIM_PTR_INT__h_s_i_s = INT__h_s_i_s + abcd * 588;
190                 double * restrict PRIM_PTR_INT__h_s_k_s = INT__h_s_k_s + abcd * 756;
191                 double * restrict PRIM_PTR_INT__i_s_g_s = INT__i_s_g_s + abcd * 420;
192                 double * restrict PRIM_PTR_INT__i_s_h_s = INT__i_s_h_s + abcd * 588;
193                 double * restrict PRIM_PTR_INT__i_s_i_s = INT__i_s_i_s + abcd * 784;
194                 double * restrict PRIM_PTR_INT__i_s_k_s = INT__i_s_k_s + abcd * 1008;
195                 double * restrict PRIM_PTR_INT__k_s_g_s = INT__k_s_g_s + abcd * 540;
196                 double * restrict PRIM_PTR_INT__k_s_h_s = INT__k_s_h_s + abcd * 756;
197                 double * restrict PRIM_PTR_INT__k_s_i_s = INT__k_s_i_s + abcd * 1008;
198                 double * restrict PRIM_PTR_INT__k_s_k_s = INT__k_s_k_s + abcd * 1296;
199 
200 
201 
202                 // Load these one per loop over i
203                 const SIMINT_DBLTYPE P_alpha = SIMINT_DBLSET1(P.alpha[i]);
204                 const SIMINT_DBLTYPE P_prefac = SIMINT_DBLSET1(P.prefac[i]);
205                 const SIMINT_DBLTYPE Pxyz[3] = { SIMINT_DBLSET1(P.x[i]), SIMINT_DBLSET1(P.y[i]), SIMINT_DBLSET1(P.z[i]) };
206 
207                 const SIMINT_DBLTYPE P_PA[3] = { SIMINT_DBLSET1(P.PA_x[i]), SIMINT_DBLSET1(P.PA_y[i]), SIMINT_DBLSET1(P.PA_z[i]) };
208 
209                 for(j = jstart; j < jend; j += SIMINT_SIMD_LEN)
210                 {
211                     // calculate the shell offsets
212                     // these are the offset from the shell pointed to by cd
213                     // for each element
214                     int shelloffsets[SIMINT_SIMD_LEN] = {0};
215                     int lastoffset = 0;
216                     const int nlane = ( ((j + SIMINT_SIMD_LEN) < jend) ? SIMINT_SIMD_LEN : (jend - j));
217 
218                     if((iprimcd + SIMINT_SIMD_LEN) >= nprim_icd)
219                     {
220                         // Handle if the first element of the vector is a new shell
221                         if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
222                         {
223                             nprim_icd += Q.nprim12[cd + (++icd)];
224                             PRIM_PTR_INT__g_s_g_s += 225;
225                             PRIM_PTR_INT__g_s_h_s += 315;
226                             PRIM_PTR_INT__g_s_i_s += 420;
227                             PRIM_PTR_INT__g_s_k_s += 540;
228                             PRIM_PTR_INT__h_s_g_s += 315;
229                             PRIM_PTR_INT__h_s_h_s += 441;
230                             PRIM_PTR_INT__h_s_i_s += 588;
231                             PRIM_PTR_INT__h_s_k_s += 756;
232                             PRIM_PTR_INT__i_s_g_s += 420;
233                             PRIM_PTR_INT__i_s_h_s += 588;
234                             PRIM_PTR_INT__i_s_i_s += 784;
235                             PRIM_PTR_INT__i_s_k_s += 1008;
236                             PRIM_PTR_INT__k_s_g_s += 540;
237                             PRIM_PTR_INT__k_s_h_s += 756;
238                             PRIM_PTR_INT__k_s_i_s += 1008;
239                             PRIM_PTR_INT__k_s_k_s += 1296;
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__g_s_g_s += lastoffset*225;
266                             PRIM_PTR_INT__g_s_h_s += lastoffset*315;
267                             PRIM_PTR_INT__g_s_i_s += lastoffset*420;
268                             PRIM_PTR_INT__g_s_k_s += lastoffset*540;
269                             PRIM_PTR_INT__h_s_g_s += lastoffset*315;
270                             PRIM_PTR_INT__h_s_h_s += lastoffset*441;
271                             PRIM_PTR_INT__h_s_i_s += lastoffset*588;
272                             PRIM_PTR_INT__h_s_k_s += lastoffset*756;
273                             PRIM_PTR_INT__i_s_g_s += lastoffset*420;
274                             PRIM_PTR_INT__i_s_h_s += lastoffset*588;
275                             PRIM_PTR_INT__i_s_i_s += lastoffset*784;
276                             PRIM_PTR_INT__i_s_k_s += lastoffset*1008;
277                             PRIM_PTR_INT__k_s_g_s += lastoffset*540;
278                             PRIM_PTR_INT__k_s_h_s += lastoffset*756;
279                             PRIM_PTR_INT__k_s_i_s += lastoffset*1008;
280                             PRIM_PTR_INT__k_s_k_s += lastoffset*1296;
281                             continue;
282                         }
283                     }
284 
285                     const SIMINT_DBLTYPE Q_alpha = SIMINT_DBLLOAD(Q.alpha, j);
286                     const SIMINT_DBLTYPE PQalpha_mul = SIMINT_MUL(P_alpha, Q_alpha);
287                     const SIMINT_DBLTYPE PQalpha_sum = SIMINT_ADD(P_alpha, Q_alpha);
288                     const SIMINT_DBLTYPE one_over_PQalpha_sum = SIMINT_DIV(const_1, PQalpha_sum);
289 
290 
291                     /* construct R2 = (Px - Qx)**2 + (Py - Qy)**2 + (Pz -Qz)**2 */
292                     SIMINT_DBLTYPE PQ[3];
293                     PQ[0] = SIMINT_SUB(Pxyz[0], SIMINT_DBLLOAD(Q.x, j));
294                     PQ[1] = SIMINT_SUB(Pxyz[1], SIMINT_DBLLOAD(Q.y, j));
295                     PQ[2] = SIMINT_SUB(Pxyz[2], SIMINT_DBLLOAD(Q.z, j));
296                     SIMINT_DBLTYPE R2 = SIMINT_MUL(PQ[0], PQ[0]);
297                     R2 = SIMINT_FMADD(PQ[1], PQ[1], R2);
298                     R2 = SIMINT_FMADD(PQ[2], PQ[2], R2);
299 
300                     const SIMINT_DBLTYPE alpha = SIMINT_MUL(PQalpha_mul, one_over_PQalpha_sum); // alpha from MEST
301                     const SIMINT_DBLTYPE one_over_p = SIMINT_DIV(const_1, P_alpha);
302                     const SIMINT_DBLTYPE one_over_q = SIMINT_DIV(const_1, Q_alpha);
303                     const SIMINT_DBLTYPE one_over_2p = SIMINT_MUL(one_half, one_over_p);
304                     const SIMINT_DBLTYPE one_over_2q = SIMINT_MUL(one_half, one_over_q);
305                     const SIMINT_DBLTYPE one_over_2pq = SIMINT_MUL(one_half, one_over_PQalpha_sum);
306                     const SIMINT_DBLTYPE Q_PA[3] = { SIMINT_DBLLOAD(Q.PA_x, j), SIMINT_DBLLOAD(Q.PA_y, j), SIMINT_DBLLOAD(Q.PA_z, j) };
307 
308                     // NOTE: Minus sign!
309                     const SIMINT_DBLTYPE a_over_p = SIMINT_MUL(SIMINT_NEG(alpha), one_over_p);
310                     SIMINT_DBLTYPE aop_PQ[3];
311                     aop_PQ[0] = SIMINT_MUL(a_over_p, PQ[0]);
312                     aop_PQ[1] = SIMINT_MUL(a_over_p, PQ[1]);
313                     aop_PQ[2] = SIMINT_MUL(a_over_p, PQ[2]);
314 
315                     SIMINT_DBLTYPE a_over_q = SIMINT_MUL(alpha, one_over_q);
316                     SIMINT_DBLTYPE aoq_PQ[3];
317                     aoq_PQ[0] = SIMINT_MUL(a_over_q, PQ[0]);
318                     aoq_PQ[1] = SIMINT_MUL(a_over_q, PQ[1]);
319                     aoq_PQ[2] = SIMINT_MUL(a_over_q, PQ[2]);
320                     // Put a minus sign here so we don't have to in RR routines
321                     a_over_q = SIMINT_NEG(a_over_q);
322 
323 
324                     //////////////////////////////////////////////
325                     // Fjt function section
326                     // Maximum v value: 14
327                     //////////////////////////////////////////////
328                     // The parameter to the Fjt function
329                     const SIMINT_DBLTYPE F_x = SIMINT_MUL(R2, alpha);
330 
331 
332                     const SIMINT_DBLTYPE Q_prefac = mask_load(nlane, Q.prefac + j);
333 
334 
335                     boys_F_split(PRIM_INT__s_s_s_s, F_x, 14);
336                     SIMINT_DBLTYPE prefac = SIMINT_SQRT(one_over_PQalpha_sum);
337                     prefac = SIMINT_MUL(SIMINT_MUL(P_prefac, Q_prefac), prefac);
338                     for(n = 0; n <= 14; n++)
339                         PRIM_INT__s_s_s_s[n] = SIMINT_MUL(PRIM_INT__s_s_s_s[n], prefac);
340 
341                     //////////////////////////////////////////////
342                     // Primitive integrals: Vertical recurrance
343                     //////////////////////////////////////////////
344 
345                     const SIMINT_DBLTYPE vrr_const_1_over_2p = one_over_2p;
346                     const SIMINT_DBLTYPE vrr_const_2_over_2p = SIMINT_MUL(const_2, one_over_2p);
347                     const SIMINT_DBLTYPE vrr_const_3_over_2p = SIMINT_MUL(const_3, one_over_2p);
348                     const SIMINT_DBLTYPE vrr_const_4_over_2p = SIMINT_MUL(const_4, one_over_2p);
349                     const SIMINT_DBLTYPE vrr_const_5_over_2p = SIMINT_MUL(const_5, one_over_2p);
350                     const SIMINT_DBLTYPE vrr_const_6_over_2p = SIMINT_MUL(const_6, one_over_2p);
351                     const SIMINT_DBLTYPE vrr_const_1_over_2q = one_over_2q;
352                     const SIMINT_DBLTYPE vrr_const_2_over_2q = SIMINT_MUL(const_2, one_over_2q);
353                     const SIMINT_DBLTYPE vrr_const_3_over_2q = SIMINT_MUL(const_3, one_over_2q);
354                     const SIMINT_DBLTYPE vrr_const_4_over_2q = SIMINT_MUL(const_4, one_over_2q);
355                     const SIMINT_DBLTYPE vrr_const_5_over_2q = SIMINT_MUL(const_5, one_over_2q);
356                     const SIMINT_DBLTYPE vrr_const_6_over_2q = SIMINT_MUL(const_6, one_over_2q);
357                     const SIMINT_DBLTYPE vrr_const_1_over_2pq = one_over_2pq;
358                     const SIMINT_DBLTYPE vrr_const_2_over_2pq = SIMINT_MUL(const_2, one_over_2pq);
359                     const SIMINT_DBLTYPE vrr_const_3_over_2pq = SIMINT_MUL(const_3, one_over_2pq);
360                     const SIMINT_DBLTYPE vrr_const_4_over_2pq = SIMINT_MUL(const_4, one_over_2pq);
361                     const SIMINT_DBLTYPE vrr_const_5_over_2pq = SIMINT_MUL(const_5, one_over_2pq);
362                     const SIMINT_DBLTYPE vrr_const_6_over_2pq = SIMINT_MUL(const_6, one_over_2pq);
363                     const SIMINT_DBLTYPE vrr_const_7_over_2pq = SIMINT_MUL(const_7, one_over_2pq);
364 
365 
366 
367                     // Forming PRIM_INT__p_s_s_s[14 * 3];
368                     for(n = 0; n < 14; ++n)  // loop over orders of auxiliary function
369                     {
370 
371                         PRIM_INT__p_s_s_s[n * 3 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
372                         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]);
373 
374                         PRIM_INT__p_s_s_s[n * 3 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
375                         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]);
376 
377                         PRIM_INT__p_s_s_s[n * 3 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
378                         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]);
379 
380                     }
381 
382 
383 
384                     // Forming PRIM_INT__d_s_s_s[13 * 6];
385                     for(n = 0; n < 13; ++n)  // loop over orders of auxiliary function
386                     {
387 
388                         PRIM_INT__d_s_s_s[n * 6 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__p_s_s_s[n * 3 + 0]);
389                         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]);
390                         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]);
391 
392                         PRIM_INT__d_s_s_s[n * 6 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__p_s_s_s[n * 3 + 0]);
393                         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]);
394 
395                         PRIM_INT__d_s_s_s[n * 6 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 0]);
396                         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]);
397 
398                         PRIM_INT__d_s_s_s[n * 6 + 3] = SIMINT_MUL(P_PA[1], PRIM_INT__p_s_s_s[n * 3 + 1]);
399                         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]);
400                         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]);
401 
402                         PRIM_INT__d_s_s_s[n * 6 + 4] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 1]);
403                         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]);
404 
405                         PRIM_INT__d_s_s_s[n * 6 + 5] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 2]);
406                         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]);
407                         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]);
408 
409                     }
410 
411 
412 
413                     // Forming PRIM_INT__f_s_s_s[12 * 10];
414                     for(n = 0; n < 12; ++n)  // loop over orders of auxiliary function
415                     {
416 
417                         PRIM_INT__f_s_s_s[n * 10 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 0]);
418                         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]);
419                         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]);
420 
421                         PRIM_INT__f_s_s_s[n * 10 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 0]);
422                         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]);
423 
424                         PRIM_INT__f_s_s_s[n * 10 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 0]);
425                         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]);
426 
427                         PRIM_INT__f_s_s_s[n * 10 + 3] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 3]);
428                         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]);
429 
430                         PRIM_INT__f_s_s_s[n * 10 + 4] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 1]);
431                         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]);
432 
433                         PRIM_INT__f_s_s_s[n * 10 + 5] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 5]);
434                         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]);
435 
436                         PRIM_INT__f_s_s_s[n * 10 + 6] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 3]);
437                         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]);
438                         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]);
439 
440                         PRIM_INT__f_s_s_s[n * 10 + 7] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 3]);
441                         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]);
442 
443                         PRIM_INT__f_s_s_s[n * 10 + 8] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 5]);
444                         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]);
445 
446                         PRIM_INT__f_s_s_s[n * 10 + 9] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 5]);
447                         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]);
448                         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]);
449 
450                     }
451 
452 
453                     VRR_I_g_s_s_s(
454                             PRIM_INT__g_s_s_s,
455                             PRIM_INT__f_s_s_s,
456                             PRIM_INT__d_s_s_s,
457                             P_PA,
458                             a_over_p,
459                             aop_PQ,
460                             one_over_2p,
461                             11);
462 
463 
464                     VRR_K_g_s_p_s(
465                             PRIM_INT__g_s_p_s,
466                             PRIM_INT__g_s_s_s,
467                             PRIM_INT__f_s_s_s,
468                             Q_PA,
469                             aoq_PQ,
470                             one_over_2pq,
471                             7);
472 
473 
474                     VRR_K_f_s_p_s(
475                             PRIM_INT__f_s_p_s,
476                             PRIM_INT__f_s_s_s,
477                             PRIM_INT__d_s_s_s,
478                             Q_PA,
479                             aoq_PQ,
480                             one_over_2pq,
481                             7);
482 
483 
484                     ostei_general_vrr_K(4, 0, 2, 0, 6,
485                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
486                             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);
487 
488 
489 
490                     // Forming PRIM_INT__d_s_p_s[7 * 18];
491                     for(n = 0; n < 7; ++n)  // loop over orders of auxiliary function
492                     {
493 
494                         PRIM_INT__d_s_p_s[n * 18 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 0]);
495                         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]);
496                         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]);
497 
498                         PRIM_INT__d_s_p_s[n * 18 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 0]);
499                         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]);
500 
501                         PRIM_INT__d_s_p_s[n * 18 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 0]);
502                         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]);
503 
504                         PRIM_INT__d_s_p_s[n * 18 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 1]);
505                         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]);
506                         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]);
507 
508                         PRIM_INT__d_s_p_s[n * 18 + 4] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 1]);
509                         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]);
510                         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]);
511 
512                         PRIM_INT__d_s_p_s[n * 18 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 1]);
513                         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]);
514 
515                         PRIM_INT__d_s_p_s[n * 18 + 6] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 2]);
516                         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]);
517                         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]);
518 
519                         PRIM_INT__d_s_p_s[n * 18 + 7] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 2]);
520                         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]);
521 
522                         PRIM_INT__d_s_p_s[n * 18 + 8] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 2]);
523                         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]);
524                         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]);
525 
526                         PRIM_INT__d_s_p_s[n * 18 + 9] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 3]);
527                         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]);
528 
529                         PRIM_INT__d_s_p_s[n * 18 + 10] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 3]);
530                         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]);
531                         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]);
532 
533                         PRIM_INT__d_s_p_s[n * 18 + 11] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 3]);
534                         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]);
535 
536                         PRIM_INT__d_s_p_s[n * 18 + 12] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 4]);
537                         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]);
538 
539                         PRIM_INT__d_s_p_s[n * 18 + 13] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 4]);
540                         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]);
541                         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]);
542 
543                         PRIM_INT__d_s_p_s[n * 18 + 14] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 4]);
544                         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]);
545                         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]);
546 
547                         PRIM_INT__d_s_p_s[n * 18 + 15] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 5]);
548                         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]);
549 
550                         PRIM_INT__d_s_p_s[n * 18 + 16] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 5]);
551                         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]);
552 
553                         PRIM_INT__d_s_p_s[n * 18 + 17] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 5]);
554                         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]);
555                         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]);
556 
557                     }
558 
559 
560                     VRR_K_f_s_d_s(
561                             PRIM_INT__f_s_d_s,
562                             PRIM_INT__f_s_p_s,
563                             PRIM_INT__f_s_s_s,
564                             PRIM_INT__d_s_p_s,
565                             Q_PA,
566                             a_over_q,
567                             aoq_PQ,
568                             one_over_2pq,
569                             one_over_2q,
570                             6);
571 
572 
573                     ostei_general_vrr_K(4, 0, 3, 0, 5,
574                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
575                             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);
576 
577 
578 
579                     // Forming PRIM_INT__p_s_p_s[7 * 9];
580                     for(n = 0; n < 7; ++n)  // loop over orders of auxiliary function
581                     {
582 
583                         PRIM_INT__p_s_p_s[n * 9 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_s_s[n * 3 + 0]);
584                         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]);
585                         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]);
586 
587                         PRIM_INT__p_s_p_s[n * 9 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_s_s[n * 3 + 0]);
588                         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]);
589 
590                         PRIM_INT__p_s_p_s[n * 9 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_s_s[n * 3 + 0]);
591                         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]);
592 
593                         PRIM_INT__p_s_p_s[n * 9 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_s_s[n * 3 + 1]);
594                         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]);
595 
596                         PRIM_INT__p_s_p_s[n * 9 + 4] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_s_s[n * 3 + 1]);
597                         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]);
598                         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]);
599 
600                         PRIM_INT__p_s_p_s[n * 9 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_s_s[n * 3 + 1]);
601                         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]);
602 
603                         PRIM_INT__p_s_p_s[n * 9 + 6] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_s_s[n * 3 + 2]);
604                         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]);
605 
606                         PRIM_INT__p_s_p_s[n * 9 + 7] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_s_s[n * 3 + 2]);
607                         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]);
608 
609                         PRIM_INT__p_s_p_s[n * 9 + 8] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_s_s[n * 3 + 2]);
610                         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]);
611                         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]);
612 
613                     }
614 
615 
616                     VRR_K_d_s_d_s(
617                             PRIM_INT__d_s_d_s,
618                             PRIM_INT__d_s_p_s,
619                             PRIM_INT__d_s_s_s,
620                             PRIM_INT__p_s_p_s,
621                             Q_PA,
622                             a_over_q,
623                             aoq_PQ,
624                             one_over_2pq,
625                             one_over_2q,
626                             6);
627 
628 
629                     ostei_general_vrr_K(3, 0, 3, 0, 5,
630                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
631                             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);
632 
633 
634                     ostei_general_vrr_K(4, 0, 4, 0, 4,
635                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
636                             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);
637 
638 
639 
640                     // Forming PRIM_INT__s_s_p_s[7 * 3];
641                     for(n = 0; n < 7; ++n)  // loop over orders of auxiliary function
642                     {
643 
644                         PRIM_INT__s_s_p_s[n * 3 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
645                         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]);
646 
647                         PRIM_INT__s_s_p_s[n * 3 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
648                         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]);
649 
650                         PRIM_INT__s_s_p_s[n * 3 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
651                         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]);
652 
653                     }
654 
655 
656 
657                     // Forming PRIM_INT__p_s_d_s[6 * 18];
658                     for(n = 0; n < 6; ++n)  // loop over orders of auxiliary function
659                     {
660 
661                         PRIM_INT__p_s_d_s[n * 18 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_p_s[n * 9 + 0]);
662                         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]);
663                         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]);
664                         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]);
665 
666                         PRIM_INT__p_s_d_s[n * 18 + 3] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_p_s[n * 9 + 1]);
667                         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]);
668                         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]);
669 
670                         PRIM_INT__p_s_d_s[n * 18 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_p_s[n * 9 + 2]);
671                         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]);
672                         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]);
673 
674                         PRIM_INT__p_s_d_s[n * 18 + 6] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_p_s[n * 9 + 3]);
675                         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]);
676                         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]);
677 
678                         PRIM_INT__p_s_d_s[n * 18 + 9] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_p_s[n * 9 + 4]);
679                         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]);
680                         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]);
681                         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]);
682 
683                         PRIM_INT__p_s_d_s[n * 18 + 11] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_p_s[n * 9 + 5]);
684                         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]);
685                         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]);
686 
687                         PRIM_INT__p_s_d_s[n * 18 + 12] = SIMINT_MUL(Q_PA[0], PRIM_INT__p_s_p_s[n * 9 + 6]);
688                         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]);
689                         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]);
690 
691                         PRIM_INT__p_s_d_s[n * 18 + 15] = SIMINT_MUL(Q_PA[1], PRIM_INT__p_s_p_s[n * 9 + 7]);
692                         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]);
693                         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]);
694 
695                         PRIM_INT__p_s_d_s[n * 18 + 17] = SIMINT_MUL(Q_PA[2], PRIM_INT__p_s_p_s[n * 9 + 8]);
696                         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]);
697                         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]);
698                         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]);
699 
700                     }
701 
702 
703                     VRR_K_d_s_f_s(
704                             PRIM_INT__d_s_f_s,
705                             PRIM_INT__d_s_d_s,
706                             PRIM_INT__d_s_p_s,
707                             PRIM_INT__p_s_d_s,
708                             Q_PA,
709                             a_over_q,
710                             aoq_PQ,
711                             one_over_2pq,
712                             one_over_2q,
713                             5);
714 
715 
716                     ostei_general_vrr_K(3, 0, 4, 0, 4,
717                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
718                             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);
719 
720 
721                     ostei_general_vrr_K(4, 0, 5, 0, 3,
722                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
723                             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);
724 
725 
726                     VRR_I_h_s_s_s(
727                             PRIM_INT__h_s_s_s,
728                             PRIM_INT__g_s_s_s,
729                             PRIM_INT__f_s_s_s,
730                             P_PA,
731                             a_over_p,
732                             aop_PQ,
733                             one_over_2p,
734                             10);
735 
736 
737                     ostei_general_vrr_K(5, 0, 1, 0, 7,
738                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
739                             PRIM_INT__h_s_s_s, NULL, NULL, PRIM_INT__g_s_s_s, NULL, PRIM_INT__h_s_p_s);
740 
741 
742                     ostei_general_vrr_K(5, 0, 2, 0, 6,
743                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
744                             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);
745 
746 
747                     ostei_general_vrr_K(5, 0, 3, 0, 5,
748                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
749                             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);
750 
751 
752                     ostei_general_vrr_K(5, 0, 4, 0, 4,
753                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
754                             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);
755 
756 
757 
758                     // Forming PRIM_INT__s_s_d_s[6 * 6];
759                     for(n = 0; n < 6; ++n)  // loop over orders of auxiliary function
760                     {
761 
762                         PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_p_s[n * 3 + 0]);
763                         PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_d_s[n * 6 + 0]);
764                         PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_s_d_s[n * 6 + 0]);
765 
766                         PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_p_s[n * 3 + 1]);
767                         PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_p_s[(n+1) * 3 + 1], PRIM_INT__s_s_d_s[n * 6 + 3]);
768                         PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_s_d_s[n * 6 + 3]);
769 
770                         PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 2]);
771                         PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_p_s[(n+1) * 3 + 2], PRIM_INT__s_s_d_s[n * 6 + 5]);
772                         PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_FMADD( vrr_const_1_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__s_s_d_s[n * 6 + 5]);
773 
774                     }
775 
776 
777                     VRR_K_p_s_f_s(
778                             PRIM_INT__p_s_f_s,
779                             PRIM_INT__p_s_d_s,
780                             PRIM_INT__p_s_p_s,
781                             PRIM_INT__s_s_d_s,
782                             Q_PA,
783                             a_over_q,
784                             aoq_PQ,
785                             one_over_2pq,
786                             one_over_2q,
787                             5);
788 
789 
790                     ostei_general_vrr_K(2, 0, 4, 0, 4,
791                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
792                             PRIM_INT__d_s_f_s, PRIM_INT__d_s_d_s, NULL, PRIM_INT__p_s_f_s, NULL, PRIM_INT__d_s_g_s);
793 
794 
795                     ostei_general_vrr_K(3, 0, 5, 0, 3,
796                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
797                             PRIM_INT__f_s_g_s, PRIM_INT__f_s_f_s, NULL, PRIM_INT__d_s_g_s, NULL, PRIM_INT__f_s_h_s);
798 
799 
800                     ostei_general_vrr_K(4, 0, 6, 0, 2,
801                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
802                             PRIM_INT__g_s_h_s, PRIM_INT__g_s_g_s, NULL, PRIM_INT__f_s_h_s, NULL, PRIM_INT__g_s_i_s);
803 
804 
805                     ostei_general_vrr_K(5, 0, 5, 0, 3,
806                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
807                             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);
808 
809 
810                     ostei_general_vrr1_I(6, 9,
811                             one_over_2p, a_over_p, aop_PQ, P_PA,
812                             PRIM_INT__h_s_s_s, PRIM_INT__g_s_s_s, PRIM_INT__i_s_s_s);
813 
814 
815                     ostei_general_vrr_K(6, 0, 1, 0, 7,
816                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
817                             PRIM_INT__i_s_s_s, NULL, NULL, PRIM_INT__h_s_s_s, NULL, PRIM_INT__i_s_p_s);
818 
819 
820                     ostei_general_vrr_K(6, 0, 2, 0, 6,
821                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
822                             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);
823 
824 
825                     ostei_general_vrr_K(6, 0, 3, 0, 5,
826                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
827                             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);
828 
829 
830                     ostei_general_vrr_K(6, 0, 4, 0, 4,
831                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
832                             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);
833 
834 
835 
836                     // Forming PRIM_INT__s_s_f_s[5 * 10];
837                     for(n = 0; n < 5; ++n)  // loop over orders of auxiliary function
838                     {
839 
840                         PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 0]);
841                         PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__s_s_f_s[n * 10 + 0]);
842                         PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_FMADD( vrr_const_2_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_p_s[n * 3 + 0]), PRIM_INT__s_s_f_s[n * 10 + 0]);
843 
844                         PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 3]);
845                         PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__s_s_f_s[n * 10 + 6]);
846                         PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_FMADD( vrr_const_2_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_p_s[(n+1) * 3 + 1], PRIM_INT__s_s_p_s[n * 3 + 1]), PRIM_INT__s_s_f_s[n * 10 + 6]);
847 
848                         PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 5]);
849                         PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__s_s_f_s[n * 10 + 9]);
850                         PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_FMADD( vrr_const_2_over_2q, SIMINT_FMADD(a_over_q, PRIM_INT__s_s_p_s[(n+1) * 3 + 2], PRIM_INT__s_s_p_s[n * 3 + 2]), PRIM_INT__s_s_f_s[n * 10 + 9]);
851 
852                     }
853 
854 
855                     VRR_K_p_s_g_s(
856                             PRIM_INT__p_s_g_s,
857                             PRIM_INT__p_s_f_s,
858                             PRIM_INT__p_s_d_s,
859                             PRIM_INT__s_s_f_s,
860                             Q_PA,
861                             a_over_q,
862                             aoq_PQ,
863                             one_over_2pq,
864                             one_over_2q,
865                             4);
866 
867 
868                     ostei_general_vrr_K(2, 0, 5, 0, 3,
869                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
870                             PRIM_INT__d_s_g_s, PRIM_INT__d_s_f_s, NULL, PRIM_INT__p_s_g_s, NULL, PRIM_INT__d_s_h_s);
871 
872 
873                     ostei_general_vrr_K(3, 0, 6, 0, 2,
874                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
875                             PRIM_INT__f_s_h_s, PRIM_INT__f_s_g_s, NULL, PRIM_INT__d_s_h_s, NULL, PRIM_INT__f_s_i_s);
876 
877 
878                     ostei_general_vrr_K(4, 0, 7, 0, 1,
879                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
880                             PRIM_INT__g_s_i_s, PRIM_INT__g_s_h_s, NULL, PRIM_INT__f_s_i_s, NULL, PRIM_INT__g_s_k_s);
881 
882 
883                     ostei_general_vrr_K(5, 0, 6, 0, 2,
884                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
885                             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);
886 
887 
888                     ostei_general_vrr_K(6, 0, 5, 0, 3,
889                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
890                             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);
891 
892 
893                     ostei_general_vrr1_I(7, 8,
894                             one_over_2p, a_over_p, aop_PQ, P_PA,
895                             PRIM_INT__i_s_s_s, PRIM_INT__h_s_s_s, PRIM_INT__k_s_s_s);
896 
897 
898                     ostei_general_vrr_K(7, 0, 1, 0, 7,
899                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
900                             PRIM_INT__k_s_s_s, NULL, NULL, PRIM_INT__i_s_s_s, NULL, PRIM_INT__k_s_p_s);
901 
902 
903                     ostei_general_vrr_K(7, 0, 2, 0, 6,
904                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
905                             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);
906 
907 
908                     ostei_general_vrr_K(7, 0, 3, 0, 5,
909                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
910                             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);
911 
912 
913                     ostei_general_vrr_K(7, 0, 4, 0, 4,
914                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
915                             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);
916 
917 
918                     ostei_general_vrr_K(5, 0, 7, 0, 1,
919                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
920                             PRIM_INT__h_s_i_s, PRIM_INT__h_s_h_s, NULL, PRIM_INT__g_s_i_s, NULL, PRIM_INT__h_s_k_s);
921 
922 
923                     ostei_general_vrr_K(6, 0, 6, 0, 2,
924                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
925                             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);
926 
927 
928                     ostei_general_vrr_K(7, 0, 5, 0, 3,
929                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
930                             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);
931 
932 
933                     ostei_general_vrr_K(6, 0, 7, 0, 1,
934                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
935                             PRIM_INT__i_s_i_s, PRIM_INT__i_s_h_s, NULL, PRIM_INT__h_s_i_s, NULL, PRIM_INT__i_s_k_s);
936 
937 
938                     ostei_general_vrr_K(7, 0, 6, 0, 2,
939                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
940                             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);
941 
942 
943                     ostei_general_vrr_K(7, 0, 7, 0, 1,
944                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
945                             PRIM_INT__k_s_i_s, PRIM_INT__k_s_h_s, NULL, PRIM_INT__i_s_i_s, NULL, PRIM_INT__k_s_k_s);
946 
947 
948 
949 
950                     ////////////////////////////////////
951                     // Accumulate contracted integrals
952                     ////////////////////////////////////
953                     if(lastoffset == 0)
954                     {
955                         contract_all(225, PRIM_INT__g_s_g_s, PRIM_PTR_INT__g_s_g_s);
956                         contract_all(315, PRIM_INT__g_s_h_s, PRIM_PTR_INT__g_s_h_s);
957                         contract_all(420, PRIM_INT__g_s_i_s, PRIM_PTR_INT__g_s_i_s);
958                         contract_all(540, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
959                         contract_all(315, PRIM_INT__h_s_g_s, PRIM_PTR_INT__h_s_g_s);
960                         contract_all(441, PRIM_INT__h_s_h_s, PRIM_PTR_INT__h_s_h_s);
961                         contract_all(588, PRIM_INT__h_s_i_s, PRIM_PTR_INT__h_s_i_s);
962                         contract_all(756, PRIM_INT__h_s_k_s, PRIM_PTR_INT__h_s_k_s);
963                         contract_all(420, PRIM_INT__i_s_g_s, PRIM_PTR_INT__i_s_g_s);
964                         contract_all(588, PRIM_INT__i_s_h_s, PRIM_PTR_INT__i_s_h_s);
965                         contract_all(784, PRIM_INT__i_s_i_s, PRIM_PTR_INT__i_s_i_s);
966                         contract_all(1008, PRIM_INT__i_s_k_s, PRIM_PTR_INT__i_s_k_s);
967                         contract_all(540, PRIM_INT__k_s_g_s, PRIM_PTR_INT__k_s_g_s);
968                         contract_all(756, PRIM_INT__k_s_h_s, PRIM_PTR_INT__k_s_h_s);
969                         contract_all(1008, PRIM_INT__k_s_i_s, PRIM_PTR_INT__k_s_i_s);
970                         contract_all(1296, PRIM_INT__k_s_k_s, PRIM_PTR_INT__k_s_k_s);
971                     }
972                     else
973                     {
974                         contract(225, shelloffsets, PRIM_INT__g_s_g_s, PRIM_PTR_INT__g_s_g_s);
975                         contract(315, shelloffsets, PRIM_INT__g_s_h_s, PRIM_PTR_INT__g_s_h_s);
976                         contract(420, shelloffsets, PRIM_INT__g_s_i_s, PRIM_PTR_INT__g_s_i_s);
977                         contract(540, shelloffsets, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
978                         contract(315, shelloffsets, PRIM_INT__h_s_g_s, PRIM_PTR_INT__h_s_g_s);
979                         contract(441, shelloffsets, PRIM_INT__h_s_h_s, PRIM_PTR_INT__h_s_h_s);
980                         contract(588, shelloffsets, PRIM_INT__h_s_i_s, PRIM_PTR_INT__h_s_i_s);
981                         contract(756, shelloffsets, PRIM_INT__h_s_k_s, PRIM_PTR_INT__h_s_k_s);
982                         contract(420, shelloffsets, PRIM_INT__i_s_g_s, PRIM_PTR_INT__i_s_g_s);
983                         contract(588, shelloffsets, PRIM_INT__i_s_h_s, PRIM_PTR_INT__i_s_h_s);
984                         contract(784, shelloffsets, PRIM_INT__i_s_i_s, PRIM_PTR_INT__i_s_i_s);
985                         contract(1008, shelloffsets, PRIM_INT__i_s_k_s, PRIM_PTR_INT__i_s_k_s);
986                         contract(540, shelloffsets, PRIM_INT__k_s_g_s, PRIM_PTR_INT__k_s_g_s);
987                         contract(756, shelloffsets, PRIM_INT__k_s_h_s, PRIM_PTR_INT__k_s_h_s);
988                         contract(1008, shelloffsets, PRIM_INT__k_s_i_s, PRIM_PTR_INT__k_s_i_s);
989                         contract(1296, shelloffsets, PRIM_INT__k_s_k_s, PRIM_PTR_INT__k_s_k_s);
990                         PRIM_PTR_INT__g_s_g_s += lastoffset*225;
991                         PRIM_PTR_INT__g_s_h_s += lastoffset*315;
992                         PRIM_PTR_INT__g_s_i_s += lastoffset*420;
993                         PRIM_PTR_INT__g_s_k_s += lastoffset*540;
994                         PRIM_PTR_INT__h_s_g_s += lastoffset*315;
995                         PRIM_PTR_INT__h_s_h_s += lastoffset*441;
996                         PRIM_PTR_INT__h_s_i_s += lastoffset*588;
997                         PRIM_PTR_INT__h_s_k_s += lastoffset*756;
998                         PRIM_PTR_INT__i_s_g_s += lastoffset*420;
999                         PRIM_PTR_INT__i_s_h_s += lastoffset*588;
1000                         PRIM_PTR_INT__i_s_i_s += lastoffset*784;
1001                         PRIM_PTR_INT__i_s_k_s += lastoffset*1008;
1002                         PRIM_PTR_INT__k_s_g_s += lastoffset*540;
1003                         PRIM_PTR_INT__k_s_h_s += lastoffset*756;
1004                         PRIM_PTR_INT__k_s_i_s += lastoffset*1008;
1005                         PRIM_PTR_INT__k_s_k_s += lastoffset*1296;
1006                     }
1007 
1008                 }  // close loop over j
1009             }  // close loop over i
1010 
1011             //Advance to the next batch
1012             jstart = SIMINT_SIMD_ROUND(jend);
1013 
1014             //////////////////////////////////////////////
1015             // Contracted integrals: Horizontal recurrance
1016             //////////////////////////////////////////////
1017 
1018 
1019             const double hAB[3] = { P.AB_x[ab], P.AB_y[ab], P.AB_z[ab] };
1020 
1021 
1022             for(abcd = 0; abcd < nshellbatch; ++abcd, ++real_abcd)
1023             {
1024                 const double hCD[3] = { Q.AB_x[cd+abcd], Q.AB_y[cd+abcd], Q.AB_z[cd+abcd] };
1025 
1026                 // set up HRR pointers
1027                 double const * restrict HRR_INT__g_s_g_s = INT__g_s_g_s + abcd * 225;
1028                 double const * restrict HRR_INT__g_s_h_s = INT__g_s_h_s + abcd * 315;
1029                 double const * restrict HRR_INT__g_s_i_s = INT__g_s_i_s + abcd * 420;
1030                 double const * restrict HRR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
1031                 double const * restrict HRR_INT__h_s_g_s = INT__h_s_g_s + abcd * 315;
1032                 double const * restrict HRR_INT__h_s_h_s = INT__h_s_h_s + abcd * 441;
1033                 double const * restrict HRR_INT__h_s_i_s = INT__h_s_i_s + abcd * 588;
1034                 double const * restrict HRR_INT__h_s_k_s = INT__h_s_k_s + abcd * 756;
1035                 double const * restrict HRR_INT__i_s_g_s = INT__i_s_g_s + abcd * 420;
1036                 double const * restrict HRR_INT__i_s_h_s = INT__i_s_h_s + abcd * 588;
1037                 double const * restrict HRR_INT__i_s_i_s = INT__i_s_i_s + abcd * 784;
1038                 double const * restrict HRR_INT__i_s_k_s = INT__i_s_k_s + abcd * 1008;
1039                 double const * restrict HRR_INT__k_s_g_s = INT__k_s_g_s + abcd * 540;
1040                 double const * restrict HRR_INT__k_s_h_s = INT__k_s_h_s + abcd * 756;
1041                 double const * restrict HRR_INT__k_s_i_s = INT__k_s_i_s + abcd * 1008;
1042                 double const * restrict HRR_INT__k_s_k_s = INT__k_s_k_s + abcd * 1296;
1043                 double * restrict HRR_INT__g_f_g_f = INT__g_f_g_f + real_abcd * 22500;
1044 
1045                 // form INT__g_p_g_s
1046                 HRR_J_g_p(
1047                     HRR_INT__g_p_g_s,
1048                     HRR_INT__g_s_g_s,
1049                     HRR_INT__h_s_g_s,
1050                     hAB, 15);
1051 
1052                 // form INT__g_p_h_s
1053                 HRR_J_g_p(
1054                     HRR_INT__g_p_h_s,
1055                     HRR_INT__g_s_h_s,
1056                     HRR_INT__h_s_h_s,
1057                     hAB, 21);
1058 
1059                 // form INT__g_p_i_s
1060                 HRR_J_g_p(
1061                     HRR_INT__g_p_i_s,
1062                     HRR_INT__g_s_i_s,
1063                     HRR_INT__h_s_i_s,
1064                     hAB, 28);
1065 
1066                 // form INT__g_p_k_s
1067                 HRR_J_g_p(
1068                     HRR_INT__g_p_k_s,
1069                     HRR_INT__g_s_k_s,
1070                     HRR_INT__h_s_k_s,
1071                     hAB, 36);
1072 
1073                 // form INT__h_p_g_s
1074                 ostei_general_hrr_J(5, 1, 4, 0, hAB, HRR_INT__i_s_g_s, HRR_INT__h_s_g_s, HRR_INT__h_p_g_s);
1075 
1076                 // form INT__h_p_h_s
1077                 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);
1078 
1079                 // form INT__h_p_i_s
1080                 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);
1081 
1082                 // form INT__h_p_k_s
1083                 ostei_general_hrr_J(5, 1, 7, 0, hAB, HRR_INT__i_s_k_s, HRR_INT__h_s_k_s, HRR_INT__h_p_k_s);
1084 
1085                 // form INT__i_p_g_s
1086                 ostei_general_hrr_J(6, 1, 4, 0, hAB, HRR_INT__k_s_g_s, HRR_INT__i_s_g_s, HRR_INT__i_p_g_s);
1087 
1088                 // form INT__i_p_h_s
1089                 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);
1090 
1091                 // form INT__i_p_i_s
1092                 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);
1093 
1094                 // form INT__i_p_k_s
1095                 ostei_general_hrr_J(6, 1, 7, 0, hAB, HRR_INT__k_s_k_s, HRR_INT__i_s_k_s, HRR_INT__i_p_k_s);
1096 
1097                 // form INT__g_d_g_s
1098                 ostei_general_hrr_J(4, 2, 4, 0, hAB, HRR_INT__h_p_g_s, HRR_INT__g_p_g_s, HRR_INT__g_d_g_s);
1099 
1100                 // form INT__g_d_h_s
1101                 ostei_general_hrr_J(4, 2, 5, 0, hAB, HRR_INT__h_p_h_s, HRR_INT__g_p_h_s, HRR_INT__g_d_h_s);
1102 
1103                 // form INT__g_d_i_s
1104                 ostei_general_hrr_J(4, 2, 6, 0, hAB, HRR_INT__h_p_i_s, HRR_INT__g_p_i_s, HRR_INT__g_d_i_s);
1105 
1106                 // form INT__g_d_k_s
1107                 ostei_general_hrr_J(4, 2, 7, 0, hAB, HRR_INT__h_p_k_s, HRR_INT__g_p_k_s, HRR_INT__g_d_k_s);
1108 
1109                 // form INT__h_d_g_s
1110                 ostei_general_hrr_J(5, 2, 4, 0, hAB, HRR_INT__i_p_g_s, HRR_INT__h_p_g_s, HRR_INT__h_d_g_s);
1111 
1112                 // form INT__h_d_h_s
1113                 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);
1114 
1115                 // form INT__h_d_i_s
1116                 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);
1117 
1118                 // form INT__h_d_k_s
1119                 ostei_general_hrr_J(5, 2, 7, 0, hAB, HRR_INT__i_p_k_s, HRR_INT__h_p_k_s, HRR_INT__h_d_k_s);
1120 
1121                 // form INT__g_f_g_s
1122                 ostei_general_hrr_J(4, 3, 4, 0, hAB, HRR_INT__h_d_g_s, HRR_INT__g_d_g_s, HRR_INT__g_f_g_s);
1123 
1124                 // form INT__g_f_h_s
1125                 ostei_general_hrr_J(4, 3, 5, 0, hAB, HRR_INT__h_d_h_s, HRR_INT__g_d_h_s, HRR_INT__g_f_h_s);
1126 
1127                 // form INT__g_f_i_s
1128                 ostei_general_hrr_J(4, 3, 6, 0, hAB, HRR_INT__h_d_i_s, HRR_INT__g_d_i_s, HRR_INT__g_f_i_s);
1129 
1130                 // form INT__g_f_k_s
1131                 ostei_general_hrr_J(4, 3, 7, 0, hAB, HRR_INT__h_d_k_s, HRR_INT__g_d_k_s, HRR_INT__g_f_k_s);
1132 
1133                 // form INT__g_f_g_p
1134                 HRR_L_g_p(
1135                     HRR_INT__g_f_g_p,
1136                     HRR_INT__g_f_g_s,
1137                     HRR_INT__g_f_h_s,
1138                     hCD, 150);
1139 
1140                 // form INT__g_f_h_p
1141                 ostei_general_hrr_L(4, 3, 5, 1, hCD, HRR_INT__g_f_i_s, HRR_INT__g_f_h_s, HRR_INT__g_f_h_p);
1142 
1143                 // form INT__g_f_i_p
1144                 ostei_general_hrr_L(4, 3, 6, 1, hCD, HRR_INT__g_f_k_s, HRR_INT__g_f_i_s, HRR_INT__g_f_i_p);
1145 
1146                 // form INT__g_f_g_d
1147                 ostei_general_hrr_L(4, 3, 4, 2, hCD, HRR_INT__g_f_h_p, HRR_INT__g_f_g_p, HRR_INT__g_f_g_d);
1148 
1149                 // form INT__g_f_h_d
1150                 ostei_general_hrr_L(4, 3, 5, 2, hCD, HRR_INT__g_f_i_p, HRR_INT__g_f_h_p, HRR_INT__g_f_h_d);
1151 
1152                 // form INT__g_f_g_f
1153                 ostei_general_hrr_L(4, 3, 4, 3, hCD, HRR_INT__g_f_h_d, HRR_INT__g_f_g_d, HRR_INT__g_f_g_f);
1154 
1155 
1156             }  // close HRR loop
1157 
1158 
1159         }   // close loop cdbatch
1160 
1161         istart = iend;
1162     }  // close loop over ab
1163 
1164     return P.nshell12_clip * Q.nshell12_clip;
1165 }
1166 
ostei_f_g_g_f(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__f_g_g_f)1167 int ostei_f_g_g_f(struct simint_multi_shellpair const P,
1168                   struct simint_multi_shellpair const Q,
1169                   double screen_tol,
1170                   double * const restrict work,
1171                   double * const restrict INT__f_g_g_f)
1172 {
1173     double P_AB[3*P.nshell12];
1174     struct simint_multi_shellpair P_tmp = P;
1175     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
1176     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
1177     P_tmp.AB_x = P_AB;
1178     P_tmp.AB_y = P_AB + P.nshell12;
1179     P_tmp.AB_z = P_AB + 2*P.nshell12;
1180 
1181     for(int i = 0; i < P.nshell12; i++)
1182     {
1183         P_tmp.AB_x[i] = -P.AB_x[i];
1184         P_tmp.AB_y[i] = -P.AB_y[i];
1185         P_tmp.AB_z[i] = -P.AB_z[i];
1186     }
1187 
1188     int ret = ostei_g_f_g_f(P_tmp, Q, screen_tol, work, INT__f_g_g_f);
1189     double buffer[22500] SIMINT_ALIGN_ARRAY_DBL;
1190 
1191     for(int q = 0; q < ret; q++)
1192     {
1193         int idx = 0;
1194         for(int a = 0; a < 10; ++a)
1195         for(int b = 0; b < 15; ++b)
1196         for(int c = 0; c < 15; ++c)
1197         for(int d = 0; d < 10; ++d)
1198             buffer[idx++] = INT__f_g_g_f[q*22500+b*1500+a*150+c*10+d];
1199 
1200         memcpy(INT__f_g_g_f+q*22500, buffer, 22500*sizeof(double));
1201     }
1202 
1203     return ret;
1204 }
1205 
ostei_g_f_f_g(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__g_f_f_g)1206 int ostei_g_f_f_g(struct simint_multi_shellpair const P,
1207                   struct simint_multi_shellpair const Q,
1208                   double screen_tol,
1209                   double * const restrict work,
1210                   double * const restrict INT__g_f_f_g)
1211 {
1212     double Q_AB[3*Q.nshell12];
1213     struct simint_multi_shellpair Q_tmp = Q;
1214     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
1215     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
1216     Q_tmp.AB_x = Q_AB;
1217     Q_tmp.AB_y = Q_AB + Q.nshell12;
1218     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
1219 
1220     for(int i = 0; i < Q.nshell12; i++)
1221     {
1222         Q_tmp.AB_x[i] = -Q.AB_x[i];
1223         Q_tmp.AB_y[i] = -Q.AB_y[i];
1224         Q_tmp.AB_z[i] = -Q.AB_z[i];
1225     }
1226 
1227     int ret = ostei_g_f_g_f(P, Q_tmp, screen_tol, work, INT__g_f_f_g);
1228     double buffer[22500] SIMINT_ALIGN_ARRAY_DBL;
1229 
1230     for(int q = 0; q < ret; q++)
1231     {
1232         int idx = 0;
1233         for(int a = 0; a < 15; ++a)
1234         for(int b = 0; b < 10; ++b)
1235         for(int c = 0; c < 10; ++c)
1236         for(int d = 0; d < 15; ++d)
1237             buffer[idx++] = INT__g_f_f_g[q*22500+a*1500+b*150+d*10+c];
1238 
1239         memcpy(INT__g_f_f_g+q*22500, buffer, 22500*sizeof(double));
1240     }
1241 
1242     return ret;
1243 }
1244 
ostei_f_g_f_g(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__f_g_f_g)1245 int ostei_f_g_f_g(struct simint_multi_shellpair const P,
1246                   struct simint_multi_shellpair const Q,
1247                   double screen_tol,
1248                   double * const restrict work,
1249                   double * const restrict INT__f_g_f_g)
1250 {
1251     double P_AB[3*P.nshell12];
1252     struct simint_multi_shellpair P_tmp = P;
1253     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
1254     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
1255     P_tmp.AB_x = P_AB;
1256     P_tmp.AB_y = P_AB + P.nshell12;
1257     P_tmp.AB_z = P_AB + 2*P.nshell12;
1258 
1259     for(int i = 0; i < P.nshell12; i++)
1260     {
1261         P_tmp.AB_x[i] = -P.AB_x[i];
1262         P_tmp.AB_y[i] = -P.AB_y[i];
1263         P_tmp.AB_z[i] = -P.AB_z[i];
1264     }
1265 
1266     double Q_AB[3*Q.nshell12];
1267     struct simint_multi_shellpair Q_tmp = Q;
1268     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
1269     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
1270     Q_tmp.AB_x = Q_AB;
1271     Q_tmp.AB_y = Q_AB + Q.nshell12;
1272     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
1273 
1274     for(int i = 0; i < Q.nshell12; i++)
1275     {
1276         Q_tmp.AB_x[i] = -Q.AB_x[i];
1277         Q_tmp.AB_y[i] = -Q.AB_y[i];
1278         Q_tmp.AB_z[i] = -Q.AB_z[i];
1279     }
1280 
1281     int ret = ostei_g_f_g_f(P_tmp, Q_tmp, screen_tol, work, INT__f_g_f_g);
1282     double buffer[22500] SIMINT_ALIGN_ARRAY_DBL;
1283 
1284     for(int q = 0; q < ret; q++)
1285     {
1286         int idx = 0;
1287         for(int a = 0; a < 10; ++a)
1288         for(int b = 0; b < 15; ++b)
1289         for(int c = 0; c < 10; ++c)
1290         for(int d = 0; d < 15; ++d)
1291             buffer[idx++] = INT__f_g_f_g[q*22500+b*1500+a*150+d*10+c];
1292 
1293         memcpy(INT__f_g_f_g+q*22500, buffer, 22500*sizeof(double));
1294     }
1295 
1296     return ret;
1297 }
1298 
1299