1 #include "simint/boys/boys.h"
2 #include "simint/ostei/gen/ostei_generated.h"
3 #include "simint/vectorization/vectorization.h"
4 #include <math.h>
5 #include <string.h>
6 
7 
ostei_f_d_k_i(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__f_d_k_i)8 int ostei_f_d_k_i(struct simint_multi_shellpair const P,
9                   struct simint_multi_shellpair const Q,
10                   double screen_tol,
11                   double * const restrict work,
12                   double * const restrict INT__f_d_k_i)
13 {
14 
15     SIMINT_ASSUME_ALIGN_DBL(work);
16     SIMINT_ASSUME_ALIGN_DBL(INT__f_d_k_i);
17     int ab, cd, abcd;
18     int istart, jstart;
19     int iprimcd, nprim_icd, icd;
20     const int check_screen = (screen_tol > 0.0);
21     int i, j;
22     int n;
23     int not_screened;
24     int real_abcd;
25     int iket;
26     int ibra;
27 
28     // partition workspace
29     double * const INT__f_s_k_s = work + (SIMINT_NSHELL_SIMD * 0);
30     double * const INT__f_s_l_s = work + (SIMINT_NSHELL_SIMD * 360);
31     double * const INT__f_s_m_s = work + (SIMINT_NSHELL_SIMD * 810);
32     double * const INT__f_s_n_s = work + (SIMINT_NSHELL_SIMD * 1360);
33     double * const INT__f_s_o_s = work + (SIMINT_NSHELL_SIMD * 2020);
34     double * const INT__f_s_q_s = work + (SIMINT_NSHELL_SIMD * 2800);
35     double * const INT__f_s_r_s = work + (SIMINT_NSHELL_SIMD * 3710);
36     double * const INT__g_s_k_s = work + (SIMINT_NSHELL_SIMD * 4760);
37     double * const INT__g_s_l_s = work + (SIMINT_NSHELL_SIMD * 5300);
38     double * const INT__g_s_m_s = work + (SIMINT_NSHELL_SIMD * 5975);
39     double * const INT__g_s_n_s = work + (SIMINT_NSHELL_SIMD * 6800);
40     double * const INT__g_s_o_s = work + (SIMINT_NSHELL_SIMD * 7790);
41     double * const INT__g_s_q_s = work + (SIMINT_NSHELL_SIMD * 8960);
42     double * const INT__g_s_r_s = work + (SIMINT_NSHELL_SIMD * 10325);
43     double * const INT__h_s_k_s = work + (SIMINT_NSHELL_SIMD * 11900);
44     double * const INT__h_s_l_s = work + (SIMINT_NSHELL_SIMD * 12656);
45     double * const INT__h_s_m_s = work + (SIMINT_NSHELL_SIMD * 13601);
46     double * const INT__h_s_n_s = work + (SIMINT_NSHELL_SIMD * 14756);
47     double * const INT__h_s_o_s = work + (SIMINT_NSHELL_SIMD * 16142);
48     double * const INT__h_s_q_s = work + (SIMINT_NSHELL_SIMD * 17780);
49     double * const INT__h_s_r_s = work + (SIMINT_NSHELL_SIMD * 19691);
50     SIMINT_DBLTYPE * const primwork = (SIMINT_DBLTYPE *)(work + SIMINT_NSHELL_SIMD*21896);
51     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_s = primwork + 0;
52     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_p_s = primwork + 19;
53     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_d_s = primwork + 73;
54     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_f_s = primwork + 175;
55     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_g_s = primwork + 335;
56     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_h_s = primwork + 560;
57     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_i_s = primwork + 854;
58     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_k_s = primwork + 1218;
59     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_l_s = primwork + 1650;
60     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_m_s = primwork + 2145;
61     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_n_s = primwork + 2695;
62     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_o_s = primwork + 3289;
63     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_q_s = primwork + 3913;
64     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_r_s = primwork + 4550;
65     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_f_s = primwork + 5180;
66     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_g_s = primwork + 5330;
67     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_h_s = primwork + 5555;
68     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_i_s = primwork + 5870;
69     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_k_s = primwork + 6290;
70     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_l_s = primwork + 6830;
71     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_m_s = primwork + 7505;
72     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_n_s = primwork + 8330;
73     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_o_s = primwork + 9320;
74     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_q_s = primwork + 10490;
75     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_r_s = primwork + 11855;
76     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_g_s = primwork + 13430;
77     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_h_s = primwork + 13790;
78     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_i_s = primwork + 14294;
79     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_k_s = primwork + 14966;
80     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_l_s = primwork + 15830;
81     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_m_s = primwork + 16910;
82     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_n_s = primwork + 18230;
83     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_o_s = primwork + 19814;
84     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_q_s = primwork + 21686;
85     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_r_s = primwork + 23870;
86     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_h_s = primwork + 26390;
87     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_i_s = primwork + 27020;
88     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_k_s = primwork + 27860;
89     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_l_s = primwork + 28940;
90     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_m_s = primwork + 30290;
91     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_n_s = primwork + 31940;
92     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_o_s = primwork + 33920;
93     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_q_s = primwork + 36260;
94     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_r_s = primwork + 38990;
95     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_i_s = primwork + 42140;
96     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_k_s = primwork + 42980;
97     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_l_s = primwork + 44060;
98     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_m_s = primwork + 45410;
99     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_n_s = primwork + 47060;
100     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_o_s = primwork + 49040;
101     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_q_s = primwork + 51380;
102     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_r_s = primwork + 54110;
103     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_k_s = primwork + 57260;
104     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_l_s = primwork + 58016;
105     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_m_s = primwork + 58961;
106     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_n_s = primwork + 60116;
107     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_o_s = primwork + 61502;
108     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_q_s = primwork + 63140;
109     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_r_s = primwork + 65051;
110     double * const hrrwork = (double *)(primwork + 67256);
111     double * const HRR_INT__f_p_k_s = hrrwork + 0;
112     double * const HRR_INT__f_p_l_s = hrrwork + 1080;
113     double * const HRR_INT__f_p_m_s = hrrwork + 2430;
114     double * const HRR_INT__f_p_n_s = hrrwork + 4080;
115     double * const HRR_INT__f_p_o_s = hrrwork + 6060;
116     double * const HRR_INT__f_p_q_s = hrrwork + 8400;
117     double * const HRR_INT__f_p_r_s = hrrwork + 11130;
118     double * const HRR_INT__f_d_k_s = hrrwork + 14280;
119     double * const HRR_INT__f_d_k_p = hrrwork + 16440;
120     double * const HRR_INT__f_d_k_d = hrrwork + 22920;
121     double * const HRR_INT__f_d_k_f = hrrwork + 35880;
122     double * const HRR_INT__f_d_k_g = hrrwork + 57480;
123     double * const HRR_INT__f_d_k_h = hrrwork + 89880;
124     double * const HRR_INT__f_d_l_s = hrrwork + 135240;
125     double * const HRR_INT__f_d_l_p = hrrwork + 137940;
126     double * const HRR_INT__f_d_l_d = hrrwork + 146040;
127     double * const HRR_INT__f_d_l_f = hrrwork + 162240;
128     double * const HRR_INT__f_d_l_g = hrrwork + 189240;
129     double * const HRR_INT__f_d_l_h = hrrwork + 229740;
130     double * const HRR_INT__f_d_m_s = hrrwork + 286440;
131     double * const HRR_INT__f_d_m_p = hrrwork + 289740;
132     double * const HRR_INT__f_d_m_d = hrrwork + 299640;
133     double * const HRR_INT__f_d_m_f = hrrwork + 319440;
134     double * const HRR_INT__f_d_m_g = hrrwork + 352440;
135     double * const HRR_INT__f_d_n_s = hrrwork + 401940;
136     double * const HRR_INT__f_d_n_p = hrrwork + 405900;
137     double * const HRR_INT__f_d_n_d = hrrwork + 417780;
138     double * const HRR_INT__f_d_n_f = hrrwork + 441540;
139     double * const HRR_INT__f_d_o_s = hrrwork + 481140;
140     double * const HRR_INT__f_d_o_p = hrrwork + 485820;
141     double * const HRR_INT__f_d_o_d = hrrwork + 499860;
142     double * const HRR_INT__f_d_q_s = hrrwork + 527940;
143     double * const HRR_INT__f_d_q_p = hrrwork + 533400;
144     double * const HRR_INT__f_d_r_s = hrrwork + 549780;
145     double * const HRR_INT__g_p_k_s = hrrwork + 556080;
146     double * const HRR_INT__g_p_l_s = hrrwork + 557700;
147     double * const HRR_INT__g_p_m_s = hrrwork + 559725;
148     double * const HRR_INT__g_p_n_s = hrrwork + 562200;
149     double * const HRR_INT__g_p_o_s = hrrwork + 565170;
150     double * const HRR_INT__g_p_q_s = hrrwork + 568680;
151     double * const HRR_INT__g_p_r_s = hrrwork + 572775;
152 
153 
154     // Create constants
155     const SIMINT_DBLTYPE const_1 = SIMINT_DBLSET1(1);
156     const SIMINT_DBLTYPE const_10 = SIMINT_DBLSET1(10);
157     const SIMINT_DBLTYPE const_11 = SIMINT_DBLSET1(11);
158     const SIMINT_DBLTYPE const_12 = SIMINT_DBLSET1(12);
159     const SIMINT_DBLTYPE const_13 = SIMINT_DBLSET1(13);
160     const SIMINT_DBLTYPE const_2 = SIMINT_DBLSET1(2);
161     const SIMINT_DBLTYPE const_3 = SIMINT_DBLSET1(3);
162     const SIMINT_DBLTYPE const_4 = SIMINT_DBLSET1(4);
163     const SIMINT_DBLTYPE const_5 = SIMINT_DBLSET1(5);
164     const SIMINT_DBLTYPE const_6 = SIMINT_DBLSET1(6);
165     const SIMINT_DBLTYPE const_7 = SIMINT_DBLSET1(7);
166     const SIMINT_DBLTYPE const_8 = SIMINT_DBLSET1(8);
167     const SIMINT_DBLTYPE const_9 = SIMINT_DBLSET1(9);
168     const SIMINT_DBLTYPE one_half = SIMINT_DBLSET1(0.5);
169 
170 
171     ////////////////////////////////////////
172     // Loop over shells and primitives
173     ////////////////////////////////////////
174 
175     real_abcd = 0;
176     istart = 0;
177     for(ab = 0; ab < P.nshell12_clip; ++ab)
178     {
179         const int iend = istart + P.nprim12[ab];
180 
181         cd = 0;
182         jstart = 0;
183 
184         for(cd = 0; cd < Q.nshell12_clip; cd += SIMINT_NSHELL_SIMD)
185         {
186             const int nshellbatch = ((cd + SIMINT_NSHELL_SIMD) > Q.nshell12_clip) ? Q.nshell12_clip - cd : SIMINT_NSHELL_SIMD;
187             int jend = jstart;
188             for(i = 0; i < nshellbatch; i++)
189                 jend += Q.nprim12[cd+i];
190 
191             // Clear the beginning of the workspace (where we are accumulating integrals)
192             memset(work, 0, SIMINT_NSHELL_SIMD * 21896 * sizeof(double));
193             abcd = 0;
194 
195 
196             for(i = istart; i < iend; ++i)
197             {
198                 SIMINT_DBLTYPE bra_screen_max;  // only used if check_screen
199 
200                 if(check_screen)
201                 {
202                     // Skip this whole thing if always insignificant
203                     if((P.screen[i] * Q.screen_max) < screen_tol)
204                         continue;
205                     bra_screen_max = SIMINT_DBLSET1(P.screen[i]);
206                 }
207 
208                 icd = 0;
209                 iprimcd = 0;
210                 nprim_icd = Q.nprim12[cd];
211                 double * restrict PRIM_PTR_INT__f_s_k_s = INT__f_s_k_s + abcd * 360;
212                 double * restrict PRIM_PTR_INT__f_s_l_s = INT__f_s_l_s + abcd * 450;
213                 double * restrict PRIM_PTR_INT__f_s_m_s = INT__f_s_m_s + abcd * 550;
214                 double * restrict PRIM_PTR_INT__f_s_n_s = INT__f_s_n_s + abcd * 660;
215                 double * restrict PRIM_PTR_INT__f_s_o_s = INT__f_s_o_s + abcd * 780;
216                 double * restrict PRIM_PTR_INT__f_s_q_s = INT__f_s_q_s + abcd * 910;
217                 double * restrict PRIM_PTR_INT__f_s_r_s = INT__f_s_r_s + abcd * 1050;
218                 double * restrict PRIM_PTR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
219                 double * restrict PRIM_PTR_INT__g_s_l_s = INT__g_s_l_s + abcd * 675;
220                 double * restrict PRIM_PTR_INT__g_s_m_s = INT__g_s_m_s + abcd * 825;
221                 double * restrict PRIM_PTR_INT__g_s_n_s = INT__g_s_n_s + abcd * 990;
222                 double * restrict PRIM_PTR_INT__g_s_o_s = INT__g_s_o_s + abcd * 1170;
223                 double * restrict PRIM_PTR_INT__g_s_q_s = INT__g_s_q_s + abcd * 1365;
224                 double * restrict PRIM_PTR_INT__g_s_r_s = INT__g_s_r_s + abcd * 1575;
225                 double * restrict PRIM_PTR_INT__h_s_k_s = INT__h_s_k_s + abcd * 756;
226                 double * restrict PRIM_PTR_INT__h_s_l_s = INT__h_s_l_s + abcd * 945;
227                 double * restrict PRIM_PTR_INT__h_s_m_s = INT__h_s_m_s + abcd * 1155;
228                 double * restrict PRIM_PTR_INT__h_s_n_s = INT__h_s_n_s + abcd * 1386;
229                 double * restrict PRIM_PTR_INT__h_s_o_s = INT__h_s_o_s + abcd * 1638;
230                 double * restrict PRIM_PTR_INT__h_s_q_s = INT__h_s_q_s + abcd * 1911;
231                 double * restrict PRIM_PTR_INT__h_s_r_s = INT__h_s_r_s + abcd * 2205;
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__f_s_k_s += 360;
258                             PRIM_PTR_INT__f_s_l_s += 450;
259                             PRIM_PTR_INT__f_s_m_s += 550;
260                             PRIM_PTR_INT__f_s_n_s += 660;
261                             PRIM_PTR_INT__f_s_o_s += 780;
262                             PRIM_PTR_INT__f_s_q_s += 910;
263                             PRIM_PTR_INT__f_s_r_s += 1050;
264                             PRIM_PTR_INT__g_s_k_s += 540;
265                             PRIM_PTR_INT__g_s_l_s += 675;
266                             PRIM_PTR_INT__g_s_m_s += 825;
267                             PRIM_PTR_INT__g_s_n_s += 990;
268                             PRIM_PTR_INT__g_s_o_s += 1170;
269                             PRIM_PTR_INT__g_s_q_s += 1365;
270                             PRIM_PTR_INT__g_s_r_s += 1575;
271                             PRIM_PTR_INT__h_s_k_s += 756;
272                             PRIM_PTR_INT__h_s_l_s += 945;
273                             PRIM_PTR_INT__h_s_m_s += 1155;
274                             PRIM_PTR_INT__h_s_n_s += 1386;
275                             PRIM_PTR_INT__h_s_o_s += 1638;
276                             PRIM_PTR_INT__h_s_q_s += 1911;
277                             PRIM_PTR_INT__h_s_r_s += 2205;
278                         }
279                         iprimcd++;
280                         for(n = 1; n < SIMINT_SIMD_LEN; ++n)
281                         {
282                             if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
283                             {
284                                 shelloffsets[n] = shelloffsets[n-1] + 1;
285                                 lastoffset++;
286                                 nprim_icd += Q.nprim12[cd + (++icd)];
287                             }
288                             else
289                                 shelloffsets[n] = shelloffsets[n-1];
290                             iprimcd++;
291                         }
292                     }
293                     else
294                         iprimcd += SIMINT_SIMD_LEN;
295 
296                     // Do we have to compute this vector (or has it been screened out)?
297                     // (not_screened != 0 means we have to do this vector)
298                     if(check_screen)
299                     {
300                         const double vmax = vector_max(SIMINT_MUL(bra_screen_max, SIMINT_DBLLOAD(Q.screen, j)));
301                         if(vmax < screen_tol)
302                         {
303                             PRIM_PTR_INT__f_s_k_s += lastoffset*360;
304                             PRIM_PTR_INT__f_s_l_s += lastoffset*450;
305                             PRIM_PTR_INT__f_s_m_s += lastoffset*550;
306                             PRIM_PTR_INT__f_s_n_s += lastoffset*660;
307                             PRIM_PTR_INT__f_s_o_s += lastoffset*780;
308                             PRIM_PTR_INT__f_s_q_s += lastoffset*910;
309                             PRIM_PTR_INT__f_s_r_s += lastoffset*1050;
310                             PRIM_PTR_INT__g_s_k_s += lastoffset*540;
311                             PRIM_PTR_INT__g_s_l_s += lastoffset*675;
312                             PRIM_PTR_INT__g_s_m_s += lastoffset*825;
313                             PRIM_PTR_INT__g_s_n_s += lastoffset*990;
314                             PRIM_PTR_INT__g_s_o_s += lastoffset*1170;
315                             PRIM_PTR_INT__g_s_q_s += lastoffset*1365;
316                             PRIM_PTR_INT__g_s_r_s += lastoffset*1575;
317                             PRIM_PTR_INT__h_s_k_s += lastoffset*756;
318                             PRIM_PTR_INT__h_s_l_s += lastoffset*945;
319                             PRIM_PTR_INT__h_s_m_s += lastoffset*1155;
320                             PRIM_PTR_INT__h_s_n_s += lastoffset*1386;
321                             PRIM_PTR_INT__h_s_o_s += lastoffset*1638;
322                             PRIM_PTR_INT__h_s_q_s += lastoffset*1911;
323                             PRIM_PTR_INT__h_s_r_s += lastoffset*2205;
324                             continue;
325                         }
326                     }
327 
328                     const SIMINT_DBLTYPE Q_alpha = SIMINT_DBLLOAD(Q.alpha, j);
329                     const SIMINT_DBLTYPE PQalpha_mul = SIMINT_MUL(P_alpha, Q_alpha);
330                     const SIMINT_DBLTYPE PQalpha_sum = SIMINT_ADD(P_alpha, Q_alpha);
331                     const SIMINT_DBLTYPE one_over_PQalpha_sum = SIMINT_DIV(const_1, PQalpha_sum);
332 
333 
334                     /* construct R2 = (Px - Qx)**2 + (Py - Qy)**2 + (Pz -Qz)**2 */
335                     SIMINT_DBLTYPE PQ[3];
336                     PQ[0] = SIMINT_SUB(Pxyz[0], SIMINT_DBLLOAD(Q.x, j));
337                     PQ[1] = SIMINT_SUB(Pxyz[1], SIMINT_DBLLOAD(Q.y, j));
338                     PQ[2] = SIMINT_SUB(Pxyz[2], SIMINT_DBLLOAD(Q.z, j));
339                     SIMINT_DBLTYPE R2 = SIMINT_MUL(PQ[0], PQ[0]);
340                     R2 = SIMINT_FMADD(PQ[1], PQ[1], R2);
341                     R2 = SIMINT_FMADD(PQ[2], PQ[2], R2);
342 
343                     const SIMINT_DBLTYPE alpha = SIMINT_MUL(PQalpha_mul, one_over_PQalpha_sum); // alpha from MEST
344                     const SIMINT_DBLTYPE one_over_p = SIMINT_DIV(const_1, P_alpha);
345                     const SIMINT_DBLTYPE one_over_q = SIMINT_DIV(const_1, Q_alpha);
346                     const SIMINT_DBLTYPE one_over_2p = SIMINT_MUL(one_half, one_over_p);
347                     const SIMINT_DBLTYPE one_over_2q = SIMINT_MUL(one_half, one_over_q);
348                     const SIMINT_DBLTYPE one_over_2pq = SIMINT_MUL(one_half, one_over_PQalpha_sum);
349                     const SIMINT_DBLTYPE Q_PA[3] = { SIMINT_DBLLOAD(Q.PA_x, j), SIMINT_DBLLOAD(Q.PA_y, j), SIMINT_DBLLOAD(Q.PA_z, j) };
350 
351                     // NOTE: Minus sign!
352                     const SIMINT_DBLTYPE a_over_p = SIMINT_MUL(SIMINT_NEG(alpha), one_over_p);
353                     SIMINT_DBLTYPE aop_PQ[3];
354                     aop_PQ[0] = SIMINT_MUL(a_over_p, PQ[0]);
355                     aop_PQ[1] = SIMINT_MUL(a_over_p, PQ[1]);
356                     aop_PQ[2] = SIMINT_MUL(a_over_p, PQ[2]);
357 
358                     SIMINT_DBLTYPE a_over_q = SIMINT_MUL(alpha, one_over_q);
359                     SIMINT_DBLTYPE aoq_PQ[3];
360                     aoq_PQ[0] = SIMINT_MUL(a_over_q, PQ[0]);
361                     aoq_PQ[1] = SIMINT_MUL(a_over_q, PQ[1]);
362                     aoq_PQ[2] = SIMINT_MUL(a_over_q, PQ[2]);
363                     // Put a minus sign here so we don't have to in RR routines
364                     a_over_q = SIMINT_NEG(a_over_q);
365 
366 
367                     //////////////////////////////////////////////
368                     // Fjt function section
369                     // Maximum v value: 18
370                     //////////////////////////////////////////////
371                     // The parameter to the Fjt function
372                     const SIMINT_DBLTYPE F_x = SIMINT_MUL(R2, alpha);
373 
374 
375                     const SIMINT_DBLTYPE Q_prefac = mask_load(nlane, Q.prefac + j);
376 
377 
378                     boys_F_split(PRIM_INT__s_s_s_s, F_x, 18);
379                     SIMINT_DBLTYPE prefac = SIMINT_SQRT(one_over_PQalpha_sum);
380                     prefac = SIMINT_MUL(SIMINT_MUL(P_prefac, Q_prefac), prefac);
381                     for(n = 0; n <= 18; n++)
382                         PRIM_INT__s_s_s_s[n] = SIMINT_MUL(PRIM_INT__s_s_s_s[n], prefac);
383 
384                     //////////////////////////////////////////////
385                     // Primitive integrals: Vertical recurrance
386                     //////////////////////////////////////////////
387 
388                     const SIMINT_DBLTYPE vrr_const_1_over_2p = one_over_2p;
389                     const SIMINT_DBLTYPE vrr_const_2_over_2p = SIMINT_MUL(const_2, one_over_2p);
390                     const SIMINT_DBLTYPE vrr_const_3_over_2p = SIMINT_MUL(const_3, one_over_2p);
391                     const SIMINT_DBLTYPE vrr_const_4_over_2p = SIMINT_MUL(const_4, one_over_2p);
392                     const SIMINT_DBLTYPE vrr_const_1_over_2q = one_over_2q;
393                     const SIMINT_DBLTYPE vrr_const_2_over_2q = SIMINT_MUL(const_2, one_over_2q);
394                     const SIMINT_DBLTYPE vrr_const_3_over_2q = SIMINT_MUL(const_3, one_over_2q);
395                     const SIMINT_DBLTYPE vrr_const_4_over_2q = SIMINT_MUL(const_4, one_over_2q);
396                     const SIMINT_DBLTYPE vrr_const_5_over_2q = SIMINT_MUL(const_5, one_over_2q);
397                     const SIMINT_DBLTYPE vrr_const_6_over_2q = SIMINT_MUL(const_6, one_over_2q);
398                     const SIMINT_DBLTYPE vrr_const_7_over_2q = SIMINT_MUL(const_7, one_over_2q);
399                     const SIMINT_DBLTYPE vrr_const_8_over_2q = SIMINT_MUL(const_8, one_over_2q);
400                     const SIMINT_DBLTYPE vrr_const_9_over_2q = SIMINT_MUL(const_9, one_over_2q);
401                     const SIMINT_DBLTYPE vrr_const_10_over_2q = SIMINT_MUL(const_10, one_over_2q);
402                     const SIMINT_DBLTYPE vrr_const_11_over_2q = SIMINT_MUL(const_11, one_over_2q);
403                     const SIMINT_DBLTYPE vrr_const_12_over_2q = SIMINT_MUL(const_12, one_over_2q);
404                     const SIMINT_DBLTYPE vrr_const_1_over_2pq = one_over_2pq;
405                     const SIMINT_DBLTYPE vrr_const_2_over_2pq = SIMINT_MUL(const_2, one_over_2pq);
406                     const SIMINT_DBLTYPE vrr_const_3_over_2pq = SIMINT_MUL(const_3, one_over_2pq);
407                     const SIMINT_DBLTYPE vrr_const_4_over_2pq = SIMINT_MUL(const_4, one_over_2pq);
408                     const SIMINT_DBLTYPE vrr_const_5_over_2pq = SIMINT_MUL(const_5, one_over_2pq);
409                     const SIMINT_DBLTYPE vrr_const_6_over_2pq = SIMINT_MUL(const_6, one_over_2pq);
410                     const SIMINT_DBLTYPE vrr_const_7_over_2pq = SIMINT_MUL(const_7, one_over_2pq);
411                     const SIMINT_DBLTYPE vrr_const_8_over_2pq = SIMINT_MUL(const_8, one_over_2pq);
412                     const SIMINT_DBLTYPE vrr_const_9_over_2pq = SIMINT_MUL(const_9, one_over_2pq);
413                     const SIMINT_DBLTYPE vrr_const_10_over_2pq = SIMINT_MUL(const_10, one_over_2pq);
414                     const SIMINT_DBLTYPE vrr_const_11_over_2pq = SIMINT_MUL(const_11, one_over_2pq);
415                     const SIMINT_DBLTYPE vrr_const_12_over_2pq = SIMINT_MUL(const_12, one_over_2pq);
416                     const SIMINT_DBLTYPE vrr_const_13_over_2pq = SIMINT_MUL(const_13, one_over_2pq);
417 
418 
419 
420                     // Forming PRIM_INT__s_s_p_s[18 * 3];
421                     for(n = 0; n < 18; ++n)  // loop over orders of auxiliary function
422                     {
423 
424                         PRIM_INT__s_s_p_s[n * 3 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
425                         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]);
426 
427                         PRIM_INT__s_s_p_s[n * 3 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
428                         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]);
429 
430                         PRIM_INT__s_s_p_s[n * 3 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
431                         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]);
432 
433                     }
434 
435 
436 
437                     // Forming PRIM_INT__s_s_d_s[17 * 6];
438                     for(n = 0; n < 17; ++n)  // loop over orders of auxiliary function
439                     {
440 
441                         PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_p_s[n * 3 + 0]);
442                         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]);
443                         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]);
444 
445                         PRIM_INT__s_s_d_s[n * 6 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_p_s[n * 3 + 0]);
446                         PRIM_INT__s_s_d_s[n * 6 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_d_s[n * 6 + 1]);
447 
448                         PRIM_INT__s_s_d_s[n * 6 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 0]);
449                         PRIM_INT__s_s_d_s[n * 6 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__s_s_d_s[n * 6 + 2]);
450 
451                         PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_p_s[n * 3 + 1]);
452                         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]);
453                         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]);
454 
455                         PRIM_INT__s_s_d_s[n * 6 + 4] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 1]);
456                         PRIM_INT__s_s_d_s[n * 6 + 4] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_p_s[(n+1) * 3 + 1], PRIM_INT__s_s_d_s[n * 6 + 4]);
457 
458                         PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 2]);
459                         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]);
460                         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]);
461 
462                     }
463 
464 
465 
466                     // Forming PRIM_INT__s_s_f_s[16 * 10];
467                     for(n = 0; n < 16; ++n)  // loop over orders of auxiliary function
468                     {
469 
470                         PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 0]);
471                         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]);
472                         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]);
473 
474                         PRIM_INT__s_s_f_s[n * 10 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 0]);
475                         PRIM_INT__s_s_f_s[n * 10 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__s_s_f_s[n * 10 + 1]);
476 
477                         PRIM_INT__s_s_f_s[n * 10 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 0]);
478                         PRIM_INT__s_s_f_s[n * 10 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__s_s_f_s[n * 10 + 2]);
479 
480                         PRIM_INT__s_s_f_s[n * 10 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 3]);
481                         PRIM_INT__s_s_f_s[n * 10 + 3] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__s_s_f_s[n * 10 + 3]);
482 
483                         PRIM_INT__s_s_f_s[n * 10 + 4] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 1]);
484                         PRIM_INT__s_s_f_s[n * 10 + 4] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 1], PRIM_INT__s_s_f_s[n * 10 + 4]);
485 
486                         PRIM_INT__s_s_f_s[n * 10 + 5] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 5]);
487                         PRIM_INT__s_s_f_s[n * 10 + 5] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__s_s_f_s[n * 10 + 5]);
488 
489                         PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 3]);
490                         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]);
491                         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]);
492 
493                         PRIM_INT__s_s_f_s[n * 10 + 7] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 3]);
494                         PRIM_INT__s_s_f_s[n * 10 + 7] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__s_s_f_s[n * 10 + 7]);
495 
496                         PRIM_INT__s_s_f_s[n * 10 + 8] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 5]);
497                         PRIM_INT__s_s_f_s[n * 10 + 8] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__s_s_f_s[n * 10 + 8]);
498 
499                         PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 5]);
500                         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]);
501                         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]);
502 
503                     }
504 
505 
506                     VRR_K_s_s_g_s(
507                             PRIM_INT__s_s_g_s,
508                             PRIM_INT__s_s_f_s,
509                             PRIM_INT__s_s_d_s,
510                             Q_PA,
511                             a_over_q,
512                             aoq_PQ,
513                             one_over_2q,
514                             15);
515 
516 
517                     VRR_K_s_s_h_s(
518                             PRIM_INT__s_s_h_s,
519                             PRIM_INT__s_s_g_s,
520                             PRIM_INT__s_s_f_s,
521                             Q_PA,
522                             a_over_q,
523                             aoq_PQ,
524                             one_over_2q,
525                             14);
526 
527 
528                     ostei_general_vrr1_K(6, 13,
529                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
530                             PRIM_INT__s_s_h_s, PRIM_INT__s_s_g_s, PRIM_INT__s_s_i_s);
531 
532 
533                     ostei_general_vrr1_K(7, 12,
534                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
535                             PRIM_INT__s_s_i_s, PRIM_INT__s_s_h_s, PRIM_INT__s_s_k_s);
536 
537 
538                     ostei_general_vrr_I(1, 0, 7, 0, 5,
539                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
540                             PRIM_INT__s_s_k_s, NULL, NULL, PRIM_INT__s_s_i_s, NULL, PRIM_INT__p_s_k_s);
541 
542 
543                     ostei_general_vrr_I(1, 0, 6, 0, 5,
544                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
545                             PRIM_INT__s_s_i_s, NULL, NULL, PRIM_INT__s_s_h_s, NULL, PRIM_INT__p_s_i_s);
546 
547 
548                     ostei_general_vrr_I(2, 0, 7, 0, 4,
549                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
550                             PRIM_INT__p_s_k_s, PRIM_INT__s_s_k_s, NULL, PRIM_INT__p_s_i_s, NULL, PRIM_INT__d_s_k_s);
551 
552 
553                     ostei_general_vrr_I(1, 0, 5, 0, 5,
554                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
555                             PRIM_INT__s_s_h_s, NULL, NULL, PRIM_INT__s_s_g_s, NULL, PRIM_INT__p_s_h_s);
556 
557 
558                     ostei_general_vrr_I(2, 0, 6, 0, 4,
559                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
560                             PRIM_INT__p_s_i_s, PRIM_INT__s_s_i_s, NULL, PRIM_INT__p_s_h_s, NULL, PRIM_INT__d_s_i_s);
561 
562 
563                     ostei_general_vrr_I(3, 0, 7, 0, 3,
564                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
565                             PRIM_INT__d_s_k_s, PRIM_INT__p_s_k_s, NULL, PRIM_INT__d_s_i_s, NULL, PRIM_INT__f_s_k_s);
566 
567 
568                     ostei_general_vrr1_K(8, 11,
569                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
570                             PRIM_INT__s_s_k_s, PRIM_INT__s_s_i_s, PRIM_INT__s_s_l_s);
571 
572 
573                     ostei_general_vrr_I(1, 0, 8, 0, 5,
574                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
575                             PRIM_INT__s_s_l_s, NULL, NULL, PRIM_INT__s_s_k_s, NULL, PRIM_INT__p_s_l_s);
576 
577 
578                     ostei_general_vrr_I(2, 0, 8, 0, 4,
579                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
580                             PRIM_INT__p_s_l_s, PRIM_INT__s_s_l_s, NULL, PRIM_INT__p_s_k_s, NULL, PRIM_INT__d_s_l_s);
581 
582 
583                     ostei_general_vrr_I(3, 0, 8, 0, 3,
584                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
585                             PRIM_INT__d_s_l_s, PRIM_INT__p_s_l_s, NULL, PRIM_INT__d_s_k_s, NULL, PRIM_INT__f_s_l_s);
586 
587 
588                     VRR_I_p_s_g_s(
589                             PRIM_INT__p_s_g_s,
590                             PRIM_INT__s_s_g_s,
591                             PRIM_INT__s_s_f_s,
592                             P_PA,
593                             aop_PQ,
594                             one_over_2pq,
595                             5);
596 
597 
598                     ostei_general_vrr_I(2, 0, 5, 0, 4,
599                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
600                             PRIM_INT__p_s_h_s, PRIM_INT__s_s_h_s, NULL, PRIM_INT__p_s_g_s, NULL, PRIM_INT__d_s_h_s);
601 
602 
603                     ostei_general_vrr_I(3, 0, 6, 0, 3,
604                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
605                             PRIM_INT__d_s_i_s, PRIM_INT__p_s_i_s, NULL, PRIM_INT__d_s_h_s, NULL, PRIM_INT__f_s_i_s);
606 
607 
608                     ostei_general_vrr_I(4, 0, 7, 0, 2,
609                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
610                             PRIM_INT__f_s_k_s, PRIM_INT__d_s_k_s, NULL, PRIM_INT__f_s_i_s, NULL, PRIM_INT__g_s_k_s);
611 
612 
613                     ostei_general_vrr1_K(9, 10,
614                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
615                             PRIM_INT__s_s_l_s, PRIM_INT__s_s_k_s, PRIM_INT__s_s_m_s);
616 
617 
618                     ostei_general_vrr_I(1, 0, 9, 0, 5,
619                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
620                             PRIM_INT__s_s_m_s, NULL, NULL, PRIM_INT__s_s_l_s, NULL, PRIM_INT__p_s_m_s);
621 
622 
623                     ostei_general_vrr_I(2, 0, 9, 0, 4,
624                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
625                             PRIM_INT__p_s_m_s, PRIM_INT__s_s_m_s, NULL, PRIM_INT__p_s_l_s, NULL, PRIM_INT__d_s_m_s);
626 
627 
628                     ostei_general_vrr_I(3, 0, 9, 0, 3,
629                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
630                             PRIM_INT__d_s_m_s, PRIM_INT__p_s_m_s, NULL, PRIM_INT__d_s_l_s, NULL, PRIM_INT__f_s_m_s);
631 
632 
633                     ostei_general_vrr_I(4, 0, 8, 0, 2,
634                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
635                             PRIM_INT__f_s_l_s, PRIM_INT__d_s_l_s, NULL, PRIM_INT__f_s_k_s, NULL, PRIM_INT__g_s_l_s);
636 
637 
638                     VRR_I_p_s_f_s(
639                             PRIM_INT__p_s_f_s,
640                             PRIM_INT__s_s_f_s,
641                             PRIM_INT__s_s_d_s,
642                             P_PA,
643                             aop_PQ,
644                             one_over_2pq,
645                             5);
646 
647 
648                     ostei_general_vrr_I(2, 0, 4, 0, 4,
649                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
650                             PRIM_INT__p_s_g_s, PRIM_INT__s_s_g_s, NULL, PRIM_INT__p_s_f_s, NULL, PRIM_INT__d_s_g_s);
651 
652 
653                     ostei_general_vrr_I(3, 0, 5, 0, 3,
654                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
655                             PRIM_INT__d_s_h_s, PRIM_INT__p_s_h_s, NULL, PRIM_INT__d_s_g_s, NULL, PRIM_INT__f_s_h_s);
656 
657 
658                     ostei_general_vrr_I(4, 0, 6, 0, 2,
659                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
660                             PRIM_INT__f_s_i_s, PRIM_INT__d_s_i_s, NULL, PRIM_INT__f_s_h_s, NULL, PRIM_INT__g_s_i_s);
661 
662 
663                     ostei_general_vrr_I(5, 0, 7, 0, 1,
664                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
665                             PRIM_INT__g_s_k_s, PRIM_INT__f_s_k_s, NULL, PRIM_INT__g_s_i_s, NULL, PRIM_INT__h_s_k_s);
666 
667 
668                     ostei_general_vrr1_K(10, 9,
669                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
670                             PRIM_INT__s_s_m_s, PRIM_INT__s_s_l_s, PRIM_INT__s_s_n_s);
671 
672 
673                     ostei_general_vrr_I(1, 0, 10, 0, 5,
674                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
675                             PRIM_INT__s_s_n_s, NULL, NULL, PRIM_INT__s_s_m_s, NULL, PRIM_INT__p_s_n_s);
676 
677 
678                     ostei_general_vrr_I(2, 0, 10, 0, 4,
679                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
680                             PRIM_INT__p_s_n_s, PRIM_INT__s_s_n_s, NULL, PRIM_INT__p_s_m_s, NULL, PRIM_INT__d_s_n_s);
681 
682 
683                     ostei_general_vrr_I(3, 0, 10, 0, 3,
684                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
685                             PRIM_INT__d_s_n_s, PRIM_INT__p_s_n_s, NULL, PRIM_INT__d_s_m_s, NULL, PRIM_INT__f_s_n_s);
686 
687 
688                     ostei_general_vrr_I(4, 0, 9, 0, 2,
689                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
690                             PRIM_INT__f_s_m_s, PRIM_INT__d_s_m_s, NULL, PRIM_INT__f_s_l_s, NULL, PRIM_INT__g_s_m_s);
691 
692 
693                     ostei_general_vrr_I(5, 0, 8, 0, 1,
694                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
695                             PRIM_INT__g_s_l_s, PRIM_INT__f_s_l_s, NULL, PRIM_INT__g_s_k_s, NULL, PRIM_INT__h_s_l_s);
696 
697 
698                     ostei_general_vrr1_K(11, 8,
699                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
700                             PRIM_INT__s_s_n_s, PRIM_INT__s_s_m_s, PRIM_INT__s_s_o_s);
701 
702 
703                     ostei_general_vrr_I(1, 0, 11, 0, 5,
704                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
705                             PRIM_INT__s_s_o_s, NULL, NULL, PRIM_INT__s_s_n_s, NULL, PRIM_INT__p_s_o_s);
706 
707 
708                     ostei_general_vrr_I(2, 0, 11, 0, 4,
709                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
710                             PRIM_INT__p_s_o_s, PRIM_INT__s_s_o_s, NULL, PRIM_INT__p_s_n_s, NULL, PRIM_INT__d_s_o_s);
711 
712 
713                     ostei_general_vrr_I(3, 0, 11, 0, 3,
714                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
715                             PRIM_INT__d_s_o_s, PRIM_INT__p_s_o_s, NULL, PRIM_INT__d_s_n_s, NULL, PRIM_INT__f_s_o_s);
716 
717 
718                     ostei_general_vrr_I(4, 0, 10, 0, 2,
719                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
720                             PRIM_INT__f_s_n_s, PRIM_INT__d_s_n_s, NULL, PRIM_INT__f_s_m_s, NULL, PRIM_INT__g_s_n_s);
721 
722 
723                     ostei_general_vrr_I(5, 0, 9, 0, 1,
724                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
725                             PRIM_INT__g_s_m_s, PRIM_INT__f_s_m_s, NULL, PRIM_INT__g_s_l_s, NULL, PRIM_INT__h_s_m_s);
726 
727 
728                     ostei_general_vrr1_K(12, 7,
729                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
730                             PRIM_INT__s_s_o_s, PRIM_INT__s_s_n_s, PRIM_INT__s_s_q_s);
731 
732 
733                     ostei_general_vrr_I(1, 0, 12, 0, 5,
734                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
735                             PRIM_INT__s_s_q_s, NULL, NULL, PRIM_INT__s_s_o_s, NULL, PRIM_INT__p_s_q_s);
736 
737 
738                     ostei_general_vrr_I(2, 0, 12, 0, 4,
739                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
740                             PRIM_INT__p_s_q_s, PRIM_INT__s_s_q_s, NULL, PRIM_INT__p_s_o_s, NULL, PRIM_INT__d_s_q_s);
741 
742 
743                     ostei_general_vrr_I(3, 0, 12, 0, 3,
744                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
745                             PRIM_INT__d_s_q_s, PRIM_INT__p_s_q_s, NULL, PRIM_INT__d_s_o_s, NULL, PRIM_INT__f_s_q_s);
746 
747 
748                     ostei_general_vrr_I(4, 0, 11, 0, 2,
749                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
750                             PRIM_INT__f_s_o_s, PRIM_INT__d_s_o_s, NULL, PRIM_INT__f_s_n_s, NULL, PRIM_INT__g_s_o_s);
751 
752 
753                     ostei_general_vrr_I(5, 0, 10, 0, 1,
754                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
755                             PRIM_INT__g_s_n_s, PRIM_INT__f_s_n_s, NULL, PRIM_INT__g_s_m_s, NULL, PRIM_INT__h_s_n_s);
756 
757 
758                     ostei_general_vrr1_K(13, 6,
759                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
760                             PRIM_INT__s_s_q_s, PRIM_INT__s_s_o_s, PRIM_INT__s_s_r_s);
761 
762 
763                     ostei_general_vrr_I(1, 0, 13, 0, 5,
764                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
765                             PRIM_INT__s_s_r_s, NULL, NULL, PRIM_INT__s_s_q_s, NULL, PRIM_INT__p_s_r_s);
766 
767 
768                     ostei_general_vrr_I(2, 0, 13, 0, 4,
769                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
770                             PRIM_INT__p_s_r_s, PRIM_INT__s_s_r_s, NULL, PRIM_INT__p_s_q_s, NULL, PRIM_INT__d_s_r_s);
771 
772 
773                     ostei_general_vrr_I(3, 0, 13, 0, 3,
774                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
775                             PRIM_INT__d_s_r_s, PRIM_INT__p_s_r_s, NULL, PRIM_INT__d_s_q_s, NULL, PRIM_INT__f_s_r_s);
776 
777 
778                     ostei_general_vrr_I(4, 0, 12, 0, 2,
779                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
780                             PRIM_INT__f_s_q_s, PRIM_INT__d_s_q_s, NULL, PRIM_INT__f_s_o_s, NULL, PRIM_INT__g_s_q_s);
781 
782 
783                     ostei_general_vrr_I(5, 0, 11, 0, 1,
784                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
785                             PRIM_INT__g_s_o_s, PRIM_INT__f_s_o_s, NULL, PRIM_INT__g_s_n_s, NULL, PRIM_INT__h_s_o_s);
786 
787 
788                     ostei_general_vrr_I(4, 0, 13, 0, 2,
789                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
790                             PRIM_INT__f_s_r_s, PRIM_INT__d_s_r_s, NULL, PRIM_INT__f_s_q_s, NULL, PRIM_INT__g_s_r_s);
791 
792 
793                     ostei_general_vrr_I(5, 0, 12, 0, 1,
794                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
795                             PRIM_INT__g_s_q_s, PRIM_INT__f_s_q_s, NULL, PRIM_INT__g_s_o_s, NULL, PRIM_INT__h_s_q_s);
796 
797 
798                     ostei_general_vrr_I(5, 0, 13, 0, 1,
799                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
800                             PRIM_INT__g_s_r_s, PRIM_INT__f_s_r_s, NULL, PRIM_INT__g_s_q_s, NULL, PRIM_INT__h_s_r_s);
801 
802 
803 
804 
805                     ////////////////////////////////////
806                     // Accumulate contracted integrals
807                     ////////////////////////////////////
808                     if(lastoffset == 0)
809                     {
810                         contract_all(360, PRIM_INT__f_s_k_s, PRIM_PTR_INT__f_s_k_s);
811                         contract_all(450, PRIM_INT__f_s_l_s, PRIM_PTR_INT__f_s_l_s);
812                         contract_all(550, PRIM_INT__f_s_m_s, PRIM_PTR_INT__f_s_m_s);
813                         contract_all(660, PRIM_INT__f_s_n_s, PRIM_PTR_INT__f_s_n_s);
814                         contract_all(780, PRIM_INT__f_s_o_s, PRIM_PTR_INT__f_s_o_s);
815                         contract_all(910, PRIM_INT__f_s_q_s, PRIM_PTR_INT__f_s_q_s);
816                         contract_all(1050, PRIM_INT__f_s_r_s, PRIM_PTR_INT__f_s_r_s);
817                         contract_all(540, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
818                         contract_all(675, PRIM_INT__g_s_l_s, PRIM_PTR_INT__g_s_l_s);
819                         contract_all(825, PRIM_INT__g_s_m_s, PRIM_PTR_INT__g_s_m_s);
820                         contract_all(990, PRIM_INT__g_s_n_s, PRIM_PTR_INT__g_s_n_s);
821                         contract_all(1170, PRIM_INT__g_s_o_s, PRIM_PTR_INT__g_s_o_s);
822                         contract_all(1365, PRIM_INT__g_s_q_s, PRIM_PTR_INT__g_s_q_s);
823                         contract_all(1575, PRIM_INT__g_s_r_s, PRIM_PTR_INT__g_s_r_s);
824                         contract_all(756, PRIM_INT__h_s_k_s, PRIM_PTR_INT__h_s_k_s);
825                         contract_all(945, PRIM_INT__h_s_l_s, PRIM_PTR_INT__h_s_l_s);
826                         contract_all(1155, PRIM_INT__h_s_m_s, PRIM_PTR_INT__h_s_m_s);
827                         contract_all(1386, PRIM_INT__h_s_n_s, PRIM_PTR_INT__h_s_n_s);
828                         contract_all(1638, PRIM_INT__h_s_o_s, PRIM_PTR_INT__h_s_o_s);
829                         contract_all(1911, PRIM_INT__h_s_q_s, PRIM_PTR_INT__h_s_q_s);
830                         contract_all(2205, PRIM_INT__h_s_r_s, PRIM_PTR_INT__h_s_r_s);
831                     }
832                     else
833                     {
834                         contract(360, shelloffsets, PRIM_INT__f_s_k_s, PRIM_PTR_INT__f_s_k_s);
835                         contract(450, shelloffsets, PRIM_INT__f_s_l_s, PRIM_PTR_INT__f_s_l_s);
836                         contract(550, shelloffsets, PRIM_INT__f_s_m_s, PRIM_PTR_INT__f_s_m_s);
837                         contract(660, shelloffsets, PRIM_INT__f_s_n_s, PRIM_PTR_INT__f_s_n_s);
838                         contract(780, shelloffsets, PRIM_INT__f_s_o_s, PRIM_PTR_INT__f_s_o_s);
839                         contract(910, shelloffsets, PRIM_INT__f_s_q_s, PRIM_PTR_INT__f_s_q_s);
840                         contract(1050, shelloffsets, PRIM_INT__f_s_r_s, PRIM_PTR_INT__f_s_r_s);
841                         contract(540, shelloffsets, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
842                         contract(675, shelloffsets, PRIM_INT__g_s_l_s, PRIM_PTR_INT__g_s_l_s);
843                         contract(825, shelloffsets, PRIM_INT__g_s_m_s, PRIM_PTR_INT__g_s_m_s);
844                         contract(990, shelloffsets, PRIM_INT__g_s_n_s, PRIM_PTR_INT__g_s_n_s);
845                         contract(1170, shelloffsets, PRIM_INT__g_s_o_s, PRIM_PTR_INT__g_s_o_s);
846                         contract(1365, shelloffsets, PRIM_INT__g_s_q_s, PRIM_PTR_INT__g_s_q_s);
847                         contract(1575, shelloffsets, PRIM_INT__g_s_r_s, PRIM_PTR_INT__g_s_r_s);
848                         contract(756, shelloffsets, PRIM_INT__h_s_k_s, PRIM_PTR_INT__h_s_k_s);
849                         contract(945, shelloffsets, PRIM_INT__h_s_l_s, PRIM_PTR_INT__h_s_l_s);
850                         contract(1155, shelloffsets, PRIM_INT__h_s_m_s, PRIM_PTR_INT__h_s_m_s);
851                         contract(1386, shelloffsets, PRIM_INT__h_s_n_s, PRIM_PTR_INT__h_s_n_s);
852                         contract(1638, shelloffsets, PRIM_INT__h_s_o_s, PRIM_PTR_INT__h_s_o_s);
853                         contract(1911, shelloffsets, PRIM_INT__h_s_q_s, PRIM_PTR_INT__h_s_q_s);
854                         contract(2205, shelloffsets, PRIM_INT__h_s_r_s, PRIM_PTR_INT__h_s_r_s);
855                         PRIM_PTR_INT__f_s_k_s += lastoffset*360;
856                         PRIM_PTR_INT__f_s_l_s += lastoffset*450;
857                         PRIM_PTR_INT__f_s_m_s += lastoffset*550;
858                         PRIM_PTR_INT__f_s_n_s += lastoffset*660;
859                         PRIM_PTR_INT__f_s_o_s += lastoffset*780;
860                         PRIM_PTR_INT__f_s_q_s += lastoffset*910;
861                         PRIM_PTR_INT__f_s_r_s += lastoffset*1050;
862                         PRIM_PTR_INT__g_s_k_s += lastoffset*540;
863                         PRIM_PTR_INT__g_s_l_s += lastoffset*675;
864                         PRIM_PTR_INT__g_s_m_s += lastoffset*825;
865                         PRIM_PTR_INT__g_s_n_s += lastoffset*990;
866                         PRIM_PTR_INT__g_s_o_s += lastoffset*1170;
867                         PRIM_PTR_INT__g_s_q_s += lastoffset*1365;
868                         PRIM_PTR_INT__g_s_r_s += lastoffset*1575;
869                         PRIM_PTR_INT__h_s_k_s += lastoffset*756;
870                         PRIM_PTR_INT__h_s_l_s += lastoffset*945;
871                         PRIM_PTR_INT__h_s_m_s += lastoffset*1155;
872                         PRIM_PTR_INT__h_s_n_s += lastoffset*1386;
873                         PRIM_PTR_INT__h_s_o_s += lastoffset*1638;
874                         PRIM_PTR_INT__h_s_q_s += lastoffset*1911;
875                         PRIM_PTR_INT__h_s_r_s += lastoffset*2205;
876                     }
877 
878                 }  // close loop over j
879             }  // close loop over i
880 
881             //Advance to the next batch
882             jstart = SIMINT_SIMD_ROUND(jend);
883 
884             //////////////////////////////////////////////
885             // Contracted integrals: Horizontal recurrance
886             //////////////////////////////////////////////
887 
888 
889             const double hAB[3] = { P.AB_x[ab], P.AB_y[ab], P.AB_z[ab] };
890 
891 
892             for(abcd = 0; abcd < nshellbatch; ++abcd, ++real_abcd)
893             {
894                 const double hCD[3] = { Q.AB_x[cd+abcd], Q.AB_y[cd+abcd], Q.AB_z[cd+abcd] };
895 
896                 // set up HRR pointers
897                 double const * restrict HRR_INT__f_s_k_s = INT__f_s_k_s + abcd * 360;
898                 double const * restrict HRR_INT__f_s_l_s = INT__f_s_l_s + abcd * 450;
899                 double const * restrict HRR_INT__f_s_m_s = INT__f_s_m_s + abcd * 550;
900                 double const * restrict HRR_INT__f_s_n_s = INT__f_s_n_s + abcd * 660;
901                 double const * restrict HRR_INT__f_s_o_s = INT__f_s_o_s + abcd * 780;
902                 double const * restrict HRR_INT__f_s_q_s = INT__f_s_q_s + abcd * 910;
903                 double const * restrict HRR_INT__f_s_r_s = INT__f_s_r_s + abcd * 1050;
904                 double const * restrict HRR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
905                 double const * restrict HRR_INT__g_s_l_s = INT__g_s_l_s + abcd * 675;
906                 double const * restrict HRR_INT__g_s_m_s = INT__g_s_m_s + abcd * 825;
907                 double const * restrict HRR_INT__g_s_n_s = INT__g_s_n_s + abcd * 990;
908                 double const * restrict HRR_INT__g_s_o_s = INT__g_s_o_s + abcd * 1170;
909                 double const * restrict HRR_INT__g_s_q_s = INT__g_s_q_s + abcd * 1365;
910                 double const * restrict HRR_INT__g_s_r_s = INT__g_s_r_s + abcd * 1575;
911                 double const * restrict HRR_INT__h_s_k_s = INT__h_s_k_s + abcd * 756;
912                 double const * restrict HRR_INT__h_s_l_s = INT__h_s_l_s + abcd * 945;
913                 double const * restrict HRR_INT__h_s_m_s = INT__h_s_m_s + abcd * 1155;
914                 double const * restrict HRR_INT__h_s_n_s = INT__h_s_n_s + abcd * 1386;
915                 double const * restrict HRR_INT__h_s_o_s = INT__h_s_o_s + abcd * 1638;
916                 double const * restrict HRR_INT__h_s_q_s = INT__h_s_q_s + abcd * 1911;
917                 double const * restrict HRR_INT__h_s_r_s = INT__h_s_r_s + abcd * 2205;
918                 double * restrict HRR_INT__f_d_k_i = INT__f_d_k_i + real_abcd * 60480;
919 
920                 // form INT__f_p_k_s
921                 HRR_J_f_p(
922                     HRR_INT__f_p_k_s,
923                     HRR_INT__f_s_k_s,
924                     HRR_INT__g_s_k_s,
925                     hAB, 36);
926 
927                 // form INT__f_p_l_s
928                 HRR_J_f_p(
929                     HRR_INT__f_p_l_s,
930                     HRR_INT__f_s_l_s,
931                     HRR_INT__g_s_l_s,
932                     hAB, 45);
933 
934                 // form INT__f_p_m_s
935                 HRR_J_f_p(
936                     HRR_INT__f_p_m_s,
937                     HRR_INT__f_s_m_s,
938                     HRR_INT__g_s_m_s,
939                     hAB, 55);
940 
941                 // form INT__f_p_n_s
942                 HRR_J_f_p(
943                     HRR_INT__f_p_n_s,
944                     HRR_INT__f_s_n_s,
945                     HRR_INT__g_s_n_s,
946                     hAB, 66);
947 
948                 // form INT__f_p_o_s
949                 HRR_J_f_p(
950                     HRR_INT__f_p_o_s,
951                     HRR_INT__f_s_o_s,
952                     HRR_INT__g_s_o_s,
953                     hAB, 78);
954 
955                 // form INT__f_p_q_s
956                 HRR_J_f_p(
957                     HRR_INT__f_p_q_s,
958                     HRR_INT__f_s_q_s,
959                     HRR_INT__g_s_q_s,
960                     hAB, 91);
961 
962                 // form INT__f_p_r_s
963                 HRR_J_f_p(
964                     HRR_INT__f_p_r_s,
965                     HRR_INT__f_s_r_s,
966                     HRR_INT__g_s_r_s,
967                     hAB, 105);
968 
969                 // form INT__g_p_k_s
970                 HRR_J_g_p(
971                     HRR_INT__g_p_k_s,
972                     HRR_INT__g_s_k_s,
973                     HRR_INT__h_s_k_s,
974                     hAB, 36);
975 
976                 // form INT__g_p_l_s
977                 HRR_J_g_p(
978                     HRR_INT__g_p_l_s,
979                     HRR_INT__g_s_l_s,
980                     HRR_INT__h_s_l_s,
981                     hAB, 45);
982 
983                 // form INT__g_p_m_s
984                 HRR_J_g_p(
985                     HRR_INT__g_p_m_s,
986                     HRR_INT__g_s_m_s,
987                     HRR_INT__h_s_m_s,
988                     hAB, 55);
989 
990                 // form INT__g_p_n_s
991                 HRR_J_g_p(
992                     HRR_INT__g_p_n_s,
993                     HRR_INT__g_s_n_s,
994                     HRR_INT__h_s_n_s,
995                     hAB, 66);
996 
997                 // form INT__g_p_o_s
998                 HRR_J_g_p(
999                     HRR_INT__g_p_o_s,
1000                     HRR_INT__g_s_o_s,
1001                     HRR_INT__h_s_o_s,
1002                     hAB, 78);
1003 
1004                 // form INT__g_p_q_s
1005                 HRR_J_g_p(
1006                     HRR_INT__g_p_q_s,
1007                     HRR_INT__g_s_q_s,
1008                     HRR_INT__h_s_q_s,
1009                     hAB, 91);
1010 
1011                 // form INT__g_p_r_s
1012                 HRR_J_g_p(
1013                     HRR_INT__g_p_r_s,
1014                     HRR_INT__g_s_r_s,
1015                     HRR_INT__h_s_r_s,
1016                     hAB, 105);
1017 
1018                 // form INT__f_d_k_s
1019                 HRR_J_f_d(
1020                     HRR_INT__f_d_k_s,
1021                     HRR_INT__f_p_k_s,
1022                     HRR_INT__g_p_k_s,
1023                     hAB, 36);
1024 
1025                 // form INT__f_d_l_s
1026                 HRR_J_f_d(
1027                     HRR_INT__f_d_l_s,
1028                     HRR_INT__f_p_l_s,
1029                     HRR_INT__g_p_l_s,
1030                     hAB, 45);
1031 
1032                 // form INT__f_d_m_s
1033                 HRR_J_f_d(
1034                     HRR_INT__f_d_m_s,
1035                     HRR_INT__f_p_m_s,
1036                     HRR_INT__g_p_m_s,
1037                     hAB, 55);
1038 
1039                 // form INT__f_d_n_s
1040                 HRR_J_f_d(
1041                     HRR_INT__f_d_n_s,
1042                     HRR_INT__f_p_n_s,
1043                     HRR_INT__g_p_n_s,
1044                     hAB, 66);
1045 
1046                 // form INT__f_d_o_s
1047                 HRR_J_f_d(
1048                     HRR_INT__f_d_o_s,
1049                     HRR_INT__f_p_o_s,
1050                     HRR_INT__g_p_o_s,
1051                     hAB, 78);
1052 
1053                 // form INT__f_d_q_s
1054                 HRR_J_f_d(
1055                     HRR_INT__f_d_q_s,
1056                     HRR_INT__f_p_q_s,
1057                     HRR_INT__g_p_q_s,
1058                     hAB, 91);
1059 
1060                 // form INT__f_d_r_s
1061                 HRR_J_f_d(
1062                     HRR_INT__f_d_r_s,
1063                     HRR_INT__f_p_r_s,
1064                     HRR_INT__g_p_r_s,
1065                     hAB, 105);
1066 
1067                 // form INT__f_d_k_p
1068                 ostei_general_hrr_L(3, 2, 7, 1, hCD, HRR_INT__f_d_l_s, HRR_INT__f_d_k_s, HRR_INT__f_d_k_p);
1069 
1070                 // form INT__f_d_l_p
1071                 ostei_general_hrr_L(3, 2, 8, 1, hCD, HRR_INT__f_d_m_s, HRR_INT__f_d_l_s, HRR_INT__f_d_l_p);
1072 
1073                 // form INT__f_d_m_p
1074                 ostei_general_hrr_L(3, 2, 9, 1, hCD, HRR_INT__f_d_n_s, HRR_INT__f_d_m_s, HRR_INT__f_d_m_p);
1075 
1076                 // form INT__f_d_n_p
1077                 ostei_general_hrr_L(3, 2, 10, 1, hCD, HRR_INT__f_d_o_s, HRR_INT__f_d_n_s, HRR_INT__f_d_n_p);
1078 
1079                 // form INT__f_d_o_p
1080                 ostei_general_hrr_L(3, 2, 11, 1, hCD, HRR_INT__f_d_q_s, HRR_INT__f_d_o_s, HRR_INT__f_d_o_p);
1081 
1082                 // form INT__f_d_q_p
1083                 ostei_general_hrr_L(3, 2, 12, 1, hCD, HRR_INT__f_d_r_s, HRR_INT__f_d_q_s, HRR_INT__f_d_q_p);
1084 
1085                 // form INT__f_d_k_d
1086                 ostei_general_hrr_L(3, 2, 7, 2, hCD, HRR_INT__f_d_l_p, HRR_INT__f_d_k_p, HRR_INT__f_d_k_d);
1087 
1088                 // form INT__f_d_l_d
1089                 ostei_general_hrr_L(3, 2, 8, 2, hCD, HRR_INT__f_d_m_p, HRR_INT__f_d_l_p, HRR_INT__f_d_l_d);
1090 
1091                 // form INT__f_d_m_d
1092                 ostei_general_hrr_L(3, 2, 9, 2, hCD, HRR_INT__f_d_n_p, HRR_INT__f_d_m_p, HRR_INT__f_d_m_d);
1093 
1094                 // form INT__f_d_n_d
1095                 ostei_general_hrr_L(3, 2, 10, 2, hCD, HRR_INT__f_d_o_p, HRR_INT__f_d_n_p, HRR_INT__f_d_n_d);
1096 
1097                 // form INT__f_d_o_d
1098                 ostei_general_hrr_L(3, 2, 11, 2, hCD, HRR_INT__f_d_q_p, HRR_INT__f_d_o_p, HRR_INT__f_d_o_d);
1099 
1100                 // form INT__f_d_k_f
1101                 ostei_general_hrr_L(3, 2, 7, 3, hCD, HRR_INT__f_d_l_d, HRR_INT__f_d_k_d, HRR_INT__f_d_k_f);
1102 
1103                 // form INT__f_d_l_f
1104                 ostei_general_hrr_L(3, 2, 8, 3, hCD, HRR_INT__f_d_m_d, HRR_INT__f_d_l_d, HRR_INT__f_d_l_f);
1105 
1106                 // form INT__f_d_m_f
1107                 ostei_general_hrr_L(3, 2, 9, 3, hCD, HRR_INT__f_d_n_d, HRR_INT__f_d_m_d, HRR_INT__f_d_m_f);
1108 
1109                 // form INT__f_d_n_f
1110                 ostei_general_hrr_L(3, 2, 10, 3, hCD, HRR_INT__f_d_o_d, HRR_INT__f_d_n_d, HRR_INT__f_d_n_f);
1111 
1112                 // form INT__f_d_k_g
1113                 ostei_general_hrr_L(3, 2, 7, 4, hCD, HRR_INT__f_d_l_f, HRR_INT__f_d_k_f, HRR_INT__f_d_k_g);
1114 
1115                 // form INT__f_d_l_g
1116                 ostei_general_hrr_L(3, 2, 8, 4, hCD, HRR_INT__f_d_m_f, HRR_INT__f_d_l_f, HRR_INT__f_d_l_g);
1117 
1118                 // form INT__f_d_m_g
1119                 ostei_general_hrr_L(3, 2, 9, 4, hCD, HRR_INT__f_d_n_f, HRR_INT__f_d_m_f, HRR_INT__f_d_m_g);
1120 
1121                 // form INT__f_d_k_h
1122                 ostei_general_hrr_L(3, 2, 7, 5, hCD, HRR_INT__f_d_l_g, HRR_INT__f_d_k_g, HRR_INT__f_d_k_h);
1123 
1124                 // form INT__f_d_l_h
1125                 ostei_general_hrr_L(3, 2, 8, 5, hCD, HRR_INT__f_d_m_g, HRR_INT__f_d_l_g, HRR_INT__f_d_l_h);
1126 
1127                 // form INT__f_d_k_i
1128                 ostei_general_hrr_L(3, 2, 7, 6, hCD, HRR_INT__f_d_l_h, HRR_INT__f_d_k_h, HRR_INT__f_d_k_i);
1129 
1130 
1131             }  // close HRR loop
1132 
1133 
1134         }   // close loop cdbatch
1135 
1136         istart = iend;
1137     }  // close loop over ab
1138 
1139     return P.nshell12_clip * Q.nshell12_clip;
1140 }
1141 
ostei_d_f_k_i(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__d_f_k_i)1142 int ostei_d_f_k_i(struct simint_multi_shellpair const P,
1143                   struct simint_multi_shellpair const Q,
1144                   double screen_tol,
1145                   double * const restrict work,
1146                   double * const restrict INT__d_f_k_i)
1147 {
1148     double P_AB[3*P.nshell12];
1149     struct simint_multi_shellpair P_tmp = P;
1150     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
1151     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
1152     P_tmp.AB_x = P_AB;
1153     P_tmp.AB_y = P_AB + P.nshell12;
1154     P_tmp.AB_z = P_AB + 2*P.nshell12;
1155 
1156     for(int i = 0; i < P.nshell12; i++)
1157     {
1158         P_tmp.AB_x[i] = -P.AB_x[i];
1159         P_tmp.AB_y[i] = -P.AB_y[i];
1160         P_tmp.AB_z[i] = -P.AB_z[i];
1161     }
1162 
1163     int ret = ostei_f_d_k_i(P_tmp, Q, screen_tol, work, INT__d_f_k_i);
1164     double buffer[60480] SIMINT_ALIGN_ARRAY_DBL;
1165 
1166     for(int q = 0; q < ret; q++)
1167     {
1168         int idx = 0;
1169         for(int a = 0; a < 6; ++a)
1170         for(int b = 0; b < 10; ++b)
1171         for(int c = 0; c < 36; ++c)
1172         for(int d = 0; d < 28; ++d)
1173             buffer[idx++] = INT__d_f_k_i[q*60480+b*6048+a*1008+c*28+d];
1174 
1175         memcpy(INT__d_f_k_i+q*60480, buffer, 60480*sizeof(double));
1176     }
1177 
1178     return ret;
1179 }
1180 
ostei_f_d_i_k(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__f_d_i_k)1181 int ostei_f_d_i_k(struct simint_multi_shellpair const P,
1182                   struct simint_multi_shellpair const Q,
1183                   double screen_tol,
1184                   double * const restrict work,
1185                   double * const restrict INT__f_d_i_k)
1186 {
1187     double Q_AB[3*Q.nshell12];
1188     struct simint_multi_shellpair Q_tmp = Q;
1189     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
1190     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
1191     Q_tmp.AB_x = Q_AB;
1192     Q_tmp.AB_y = Q_AB + Q.nshell12;
1193     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
1194 
1195     for(int i = 0; i < Q.nshell12; i++)
1196     {
1197         Q_tmp.AB_x[i] = -Q.AB_x[i];
1198         Q_tmp.AB_y[i] = -Q.AB_y[i];
1199         Q_tmp.AB_z[i] = -Q.AB_z[i];
1200     }
1201 
1202     int ret = ostei_f_d_k_i(P, Q_tmp, screen_tol, work, INT__f_d_i_k);
1203     double buffer[60480] SIMINT_ALIGN_ARRAY_DBL;
1204 
1205     for(int q = 0; q < ret; q++)
1206     {
1207         int idx = 0;
1208         for(int a = 0; a < 10; ++a)
1209         for(int b = 0; b < 6; ++b)
1210         for(int c = 0; c < 28; ++c)
1211         for(int d = 0; d < 36; ++d)
1212             buffer[idx++] = INT__f_d_i_k[q*60480+a*6048+b*1008+d*28+c];
1213 
1214         memcpy(INT__f_d_i_k+q*60480, buffer, 60480*sizeof(double));
1215     }
1216 
1217     return ret;
1218 }
1219 
ostei_d_f_i_k(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__d_f_i_k)1220 int ostei_d_f_i_k(struct simint_multi_shellpair const P,
1221                   struct simint_multi_shellpair const Q,
1222                   double screen_tol,
1223                   double * const restrict work,
1224                   double * const restrict INT__d_f_i_k)
1225 {
1226     double P_AB[3*P.nshell12];
1227     struct simint_multi_shellpair P_tmp = P;
1228     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
1229     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
1230     P_tmp.AB_x = P_AB;
1231     P_tmp.AB_y = P_AB + P.nshell12;
1232     P_tmp.AB_z = P_AB + 2*P.nshell12;
1233 
1234     for(int i = 0; i < P.nshell12; i++)
1235     {
1236         P_tmp.AB_x[i] = -P.AB_x[i];
1237         P_tmp.AB_y[i] = -P.AB_y[i];
1238         P_tmp.AB_z[i] = -P.AB_z[i];
1239     }
1240 
1241     double Q_AB[3*Q.nshell12];
1242     struct simint_multi_shellpair Q_tmp = Q;
1243     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
1244     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
1245     Q_tmp.AB_x = Q_AB;
1246     Q_tmp.AB_y = Q_AB + Q.nshell12;
1247     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
1248 
1249     for(int i = 0; i < Q.nshell12; i++)
1250     {
1251         Q_tmp.AB_x[i] = -Q.AB_x[i];
1252         Q_tmp.AB_y[i] = -Q.AB_y[i];
1253         Q_tmp.AB_z[i] = -Q.AB_z[i];
1254     }
1255 
1256     int ret = ostei_f_d_k_i(P_tmp, Q_tmp, screen_tol, work, INT__d_f_i_k);
1257     double buffer[60480] SIMINT_ALIGN_ARRAY_DBL;
1258 
1259     for(int q = 0; q < ret; q++)
1260     {
1261         int idx = 0;
1262         for(int a = 0; a < 6; ++a)
1263         for(int b = 0; b < 10; ++b)
1264         for(int c = 0; c < 28; ++c)
1265         for(int d = 0; d < 36; ++d)
1266             buffer[idx++] = INT__d_f_i_k[q*60480+b*6048+a*1008+d*28+c];
1267 
1268         memcpy(INT__d_f_i_k+q*60480, buffer, 60480*sizeof(double));
1269     }
1270 
1271     return ret;
1272 }
1273 
1274