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