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_p_i_h(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__g_p_i_h)8 int ostei_g_p_i_h(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_p_i_h)
13 {
14 
15     SIMINT_ASSUME_ALIGN_DBL(work);
16     SIMINT_ASSUME_ALIGN_DBL(INT__g_p_i_h);
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_i_s = work + (SIMINT_NSHELL_SIMD * 0);
30     double * const INT__g_s_k_s = work + (SIMINT_NSHELL_SIMD * 420);
31     double * const INT__g_s_l_s = work + (SIMINT_NSHELL_SIMD * 960);
32     double * const INT__g_s_m_s = work + (SIMINT_NSHELL_SIMD * 1635);
33     double * const INT__g_s_n_s = work + (SIMINT_NSHELL_SIMD * 2460);
34     double * const INT__g_s_o_s = work + (SIMINT_NSHELL_SIMD * 3450);
35     double * const INT__h_s_i_s = work + (SIMINT_NSHELL_SIMD * 4620);
36     double * const INT__h_s_k_s = work + (SIMINT_NSHELL_SIMD * 5208);
37     double * const INT__h_s_l_s = work + (SIMINT_NSHELL_SIMD * 5964);
38     double * const INT__h_s_m_s = work + (SIMINT_NSHELL_SIMD * 6909);
39     double * const INT__h_s_n_s = work + (SIMINT_NSHELL_SIMD * 8064);
40     double * const INT__h_s_o_s = work + (SIMINT_NSHELL_SIMD * 9450);
41     SIMINT_DBLTYPE * const primwork = (SIMINT_DBLTYPE *)(work + SIMINT_NSHELL_SIMD*11088);
42     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_s = primwork + 0;
43     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_p_s = primwork + 17;
44     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_d_s = primwork + 65;
45     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_f_s = primwork + 155;
46     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_g_s = primwork + 295;
47     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_h_s = primwork + 490;
48     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_i_s = primwork + 742;
49     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_k_s = primwork + 1050;
50     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_l_s = primwork + 1410;
51     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_m_s = primwork + 1815;
52     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_n_s = primwork + 2255;
53     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_o_s = primwork + 2717;
54     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_d_s = primwork + 3185;
55     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_f_s = primwork + 3275;
56     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_g_s = primwork + 3425;
57     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_h_s = primwork + 3650;
58     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_i_s = primwork + 3965;
59     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_k_s = primwork + 4385;
60     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_l_s = primwork + 4925;
61     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_m_s = primwork + 5600;
62     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_n_s = primwork + 6425;
63     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_o_s = primwork + 7415;
64     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_f_s = primwork + 8585;
65     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_g_s = primwork + 8825;
66     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_h_s = primwork + 9185;
67     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_i_s = primwork + 9689;
68     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_k_s = primwork + 10361;
69     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_l_s = primwork + 11225;
70     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_m_s = primwork + 12305;
71     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_n_s = primwork + 13625;
72     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_o_s = primwork + 15209;
73     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_g_s = primwork + 17081;
74     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_h_s = primwork + 17531;
75     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_i_s = primwork + 18161;
76     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_k_s = primwork + 19001;
77     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_l_s = primwork + 20081;
78     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_m_s = primwork + 21431;
79     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_n_s = primwork + 23081;
80     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_o_s = primwork + 25061;
81     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_h_s = primwork + 27401;
82     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_i_s = primwork + 28031;
83     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_k_s = primwork + 28871;
84     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_l_s = primwork + 29951;
85     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_m_s = primwork + 31301;
86     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_n_s = primwork + 32951;
87     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_o_s = primwork + 34931;
88     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_i_s = primwork + 37271;
89     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_k_s = primwork + 37859;
90     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_l_s = primwork + 38615;
91     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_m_s = primwork + 39560;
92     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_n_s = primwork + 40715;
93     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_o_s = primwork + 42101;
94     double * const hrrwork = (double *)(primwork + 43739);
95     double * const HRR_INT__g_p_i_s = hrrwork + 0;
96     double * const HRR_INT__g_p_i_p = hrrwork + 1260;
97     double * const HRR_INT__g_p_i_d = hrrwork + 5040;
98     double * const HRR_INT__g_p_i_f = hrrwork + 12600;
99     double * const HRR_INT__g_p_i_g = hrrwork + 25200;
100     double * const HRR_INT__g_p_k_s = hrrwork + 44100;
101     double * const HRR_INT__g_p_k_p = hrrwork + 45720;
102     double * const HRR_INT__g_p_k_d = hrrwork + 50580;
103     double * const HRR_INT__g_p_k_f = hrrwork + 60300;
104     double * const HRR_INT__g_p_k_g = hrrwork + 76500;
105     double * const HRR_INT__g_p_l_s = hrrwork + 100800;
106     double * const HRR_INT__g_p_l_p = hrrwork + 102825;
107     double * const HRR_INT__g_p_l_d = hrrwork + 108900;
108     double * const HRR_INT__g_p_l_f = hrrwork + 121050;
109     double * const HRR_INT__g_p_m_s = hrrwork + 141300;
110     double * const HRR_INT__g_p_m_p = hrrwork + 143775;
111     double * const HRR_INT__g_p_m_d = hrrwork + 151200;
112     double * const HRR_INT__g_p_n_s = hrrwork + 166050;
113     double * const HRR_INT__g_p_n_p = hrrwork + 169020;
114     double * const HRR_INT__g_p_o_s = hrrwork + 177930;
115 
116 
117     // Create constants
118     const SIMINT_DBLTYPE const_1 = SIMINT_DBLSET1(1);
119     const SIMINT_DBLTYPE const_10 = SIMINT_DBLSET1(10);
120     const SIMINT_DBLTYPE const_11 = SIMINT_DBLSET1(11);
121     const SIMINT_DBLTYPE const_2 = SIMINT_DBLSET1(2);
122     const SIMINT_DBLTYPE const_3 = SIMINT_DBLSET1(3);
123     const SIMINT_DBLTYPE const_4 = SIMINT_DBLSET1(4);
124     const SIMINT_DBLTYPE const_5 = SIMINT_DBLSET1(5);
125     const SIMINT_DBLTYPE const_6 = SIMINT_DBLSET1(6);
126     const SIMINT_DBLTYPE const_7 = SIMINT_DBLSET1(7);
127     const SIMINT_DBLTYPE const_8 = SIMINT_DBLSET1(8);
128     const SIMINT_DBLTYPE const_9 = SIMINT_DBLSET1(9);
129     const SIMINT_DBLTYPE one_half = SIMINT_DBLSET1(0.5);
130 
131 
132     ////////////////////////////////////////
133     // Loop over shells and primitives
134     ////////////////////////////////////////
135 
136     real_abcd = 0;
137     istart = 0;
138     for(ab = 0; ab < P.nshell12_clip; ++ab)
139     {
140         const int iend = istart + P.nprim12[ab];
141 
142         cd = 0;
143         jstart = 0;
144 
145         for(cd = 0; cd < Q.nshell12_clip; cd += SIMINT_NSHELL_SIMD)
146         {
147             const int nshellbatch = ((cd + SIMINT_NSHELL_SIMD) > Q.nshell12_clip) ? Q.nshell12_clip - cd : SIMINT_NSHELL_SIMD;
148             int jend = jstart;
149             for(i = 0; i < nshellbatch; i++)
150                 jend += Q.nprim12[cd+i];
151 
152             // Clear the beginning of the workspace (where we are accumulating integrals)
153             memset(work, 0, SIMINT_NSHELL_SIMD * 11088 * sizeof(double));
154             abcd = 0;
155 
156 
157             for(i = istart; i < iend; ++i)
158             {
159                 SIMINT_DBLTYPE bra_screen_max;  // only used if check_screen
160 
161                 if(check_screen)
162                 {
163                     // Skip this whole thing if always insignificant
164                     if((P.screen[i] * Q.screen_max) < screen_tol)
165                         continue;
166                     bra_screen_max = SIMINT_DBLSET1(P.screen[i]);
167                 }
168 
169                 icd = 0;
170                 iprimcd = 0;
171                 nprim_icd = Q.nprim12[cd];
172                 double * restrict PRIM_PTR_INT__g_s_i_s = INT__g_s_i_s + abcd * 420;
173                 double * restrict PRIM_PTR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
174                 double * restrict PRIM_PTR_INT__g_s_l_s = INT__g_s_l_s + abcd * 675;
175                 double * restrict PRIM_PTR_INT__g_s_m_s = INT__g_s_m_s + abcd * 825;
176                 double * restrict PRIM_PTR_INT__g_s_n_s = INT__g_s_n_s + abcd * 990;
177                 double * restrict PRIM_PTR_INT__g_s_o_s = INT__g_s_o_s + abcd * 1170;
178                 double * restrict PRIM_PTR_INT__h_s_i_s = INT__h_s_i_s + abcd * 588;
179                 double * restrict PRIM_PTR_INT__h_s_k_s = INT__h_s_k_s + abcd * 756;
180                 double * restrict PRIM_PTR_INT__h_s_l_s = INT__h_s_l_s + abcd * 945;
181                 double * restrict PRIM_PTR_INT__h_s_m_s = INT__h_s_m_s + abcd * 1155;
182                 double * restrict PRIM_PTR_INT__h_s_n_s = INT__h_s_n_s + abcd * 1386;
183                 double * restrict PRIM_PTR_INT__h_s_o_s = INT__h_s_o_s + abcd * 1638;
184 
185 
186 
187                 // Load these one per loop over i
188                 const SIMINT_DBLTYPE P_alpha = SIMINT_DBLSET1(P.alpha[i]);
189                 const SIMINT_DBLTYPE P_prefac = SIMINT_DBLSET1(P.prefac[i]);
190                 const SIMINT_DBLTYPE Pxyz[3] = { SIMINT_DBLSET1(P.x[i]), SIMINT_DBLSET1(P.y[i]), SIMINT_DBLSET1(P.z[i]) };
191 
192                 const SIMINT_DBLTYPE P_PA[3] = { SIMINT_DBLSET1(P.PA_x[i]), SIMINT_DBLSET1(P.PA_y[i]), SIMINT_DBLSET1(P.PA_z[i]) };
193 
194                 for(j = jstart; j < jend; j += SIMINT_SIMD_LEN)
195                 {
196                     // calculate the shell offsets
197                     // these are the offset from the shell pointed to by cd
198                     // for each element
199                     int shelloffsets[SIMINT_SIMD_LEN] = {0};
200                     int lastoffset = 0;
201                     const int nlane = ( ((j + SIMINT_SIMD_LEN) < jend) ? SIMINT_SIMD_LEN : (jend - j));
202 
203                     if((iprimcd + SIMINT_SIMD_LEN) >= nprim_icd)
204                     {
205                         // Handle if the first element of the vector is a new shell
206                         if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
207                         {
208                             nprim_icd += Q.nprim12[cd + (++icd)];
209                             PRIM_PTR_INT__g_s_i_s += 420;
210                             PRIM_PTR_INT__g_s_k_s += 540;
211                             PRIM_PTR_INT__g_s_l_s += 675;
212                             PRIM_PTR_INT__g_s_m_s += 825;
213                             PRIM_PTR_INT__g_s_n_s += 990;
214                             PRIM_PTR_INT__g_s_o_s += 1170;
215                             PRIM_PTR_INT__h_s_i_s += 588;
216                             PRIM_PTR_INT__h_s_k_s += 756;
217                             PRIM_PTR_INT__h_s_l_s += 945;
218                             PRIM_PTR_INT__h_s_m_s += 1155;
219                             PRIM_PTR_INT__h_s_n_s += 1386;
220                             PRIM_PTR_INT__h_s_o_s += 1638;
221                         }
222                         iprimcd++;
223                         for(n = 1; n < SIMINT_SIMD_LEN; ++n)
224                         {
225                             if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
226                             {
227                                 shelloffsets[n] = shelloffsets[n-1] + 1;
228                                 lastoffset++;
229                                 nprim_icd += Q.nprim12[cd + (++icd)];
230                             }
231                             else
232                                 shelloffsets[n] = shelloffsets[n-1];
233                             iprimcd++;
234                         }
235                     }
236                     else
237                         iprimcd += SIMINT_SIMD_LEN;
238 
239                     // Do we have to compute this vector (or has it been screened out)?
240                     // (not_screened != 0 means we have to do this vector)
241                     if(check_screen)
242                     {
243                         const double vmax = vector_max(SIMINT_MUL(bra_screen_max, SIMINT_DBLLOAD(Q.screen, j)));
244                         if(vmax < screen_tol)
245                         {
246                             PRIM_PTR_INT__g_s_i_s += lastoffset*420;
247                             PRIM_PTR_INT__g_s_k_s += lastoffset*540;
248                             PRIM_PTR_INT__g_s_l_s += lastoffset*675;
249                             PRIM_PTR_INT__g_s_m_s += lastoffset*825;
250                             PRIM_PTR_INT__g_s_n_s += lastoffset*990;
251                             PRIM_PTR_INT__g_s_o_s += lastoffset*1170;
252                             PRIM_PTR_INT__h_s_i_s += lastoffset*588;
253                             PRIM_PTR_INT__h_s_k_s += lastoffset*756;
254                             PRIM_PTR_INT__h_s_l_s += lastoffset*945;
255                             PRIM_PTR_INT__h_s_m_s += lastoffset*1155;
256                             PRIM_PTR_INT__h_s_n_s += lastoffset*1386;
257                             PRIM_PTR_INT__h_s_o_s += lastoffset*1638;
258                             continue;
259                         }
260                     }
261 
262                     const SIMINT_DBLTYPE Q_alpha = SIMINT_DBLLOAD(Q.alpha, j);
263                     const SIMINT_DBLTYPE PQalpha_mul = SIMINT_MUL(P_alpha, Q_alpha);
264                     const SIMINT_DBLTYPE PQalpha_sum = SIMINT_ADD(P_alpha, Q_alpha);
265                     const SIMINT_DBLTYPE one_over_PQalpha_sum = SIMINT_DIV(const_1, PQalpha_sum);
266 
267 
268                     /* construct R2 = (Px - Qx)**2 + (Py - Qy)**2 + (Pz -Qz)**2 */
269                     SIMINT_DBLTYPE PQ[3];
270                     PQ[0] = SIMINT_SUB(Pxyz[0], SIMINT_DBLLOAD(Q.x, j));
271                     PQ[1] = SIMINT_SUB(Pxyz[1], SIMINT_DBLLOAD(Q.y, j));
272                     PQ[2] = SIMINT_SUB(Pxyz[2], SIMINT_DBLLOAD(Q.z, j));
273                     SIMINT_DBLTYPE R2 = SIMINT_MUL(PQ[0], PQ[0]);
274                     R2 = SIMINT_FMADD(PQ[1], PQ[1], R2);
275                     R2 = SIMINT_FMADD(PQ[2], PQ[2], R2);
276 
277                     const SIMINT_DBLTYPE alpha = SIMINT_MUL(PQalpha_mul, one_over_PQalpha_sum); // alpha from MEST
278                     const SIMINT_DBLTYPE one_over_p = SIMINT_DIV(const_1, P_alpha);
279                     const SIMINT_DBLTYPE one_over_q = SIMINT_DIV(const_1, Q_alpha);
280                     const SIMINT_DBLTYPE one_over_2p = SIMINT_MUL(one_half, one_over_p);
281                     const SIMINT_DBLTYPE one_over_2q = SIMINT_MUL(one_half, one_over_q);
282                     const SIMINT_DBLTYPE one_over_2pq = SIMINT_MUL(one_half, one_over_PQalpha_sum);
283                     const SIMINT_DBLTYPE Q_PA[3] = { SIMINT_DBLLOAD(Q.PA_x, j), SIMINT_DBLLOAD(Q.PA_y, j), SIMINT_DBLLOAD(Q.PA_z, j) };
284 
285                     // NOTE: Minus sign!
286                     const SIMINT_DBLTYPE a_over_p = SIMINT_MUL(SIMINT_NEG(alpha), one_over_p);
287                     SIMINT_DBLTYPE aop_PQ[3];
288                     aop_PQ[0] = SIMINT_MUL(a_over_p, PQ[0]);
289                     aop_PQ[1] = SIMINT_MUL(a_over_p, PQ[1]);
290                     aop_PQ[2] = SIMINT_MUL(a_over_p, PQ[2]);
291 
292                     SIMINT_DBLTYPE a_over_q = SIMINT_MUL(alpha, one_over_q);
293                     SIMINT_DBLTYPE aoq_PQ[3];
294                     aoq_PQ[0] = SIMINT_MUL(a_over_q, PQ[0]);
295                     aoq_PQ[1] = SIMINT_MUL(a_over_q, PQ[1]);
296                     aoq_PQ[2] = SIMINT_MUL(a_over_q, PQ[2]);
297                     // Put a minus sign here so we don't have to in RR routines
298                     a_over_q = SIMINT_NEG(a_over_q);
299 
300 
301                     //////////////////////////////////////////////
302                     // Fjt function section
303                     // Maximum v value: 16
304                     //////////////////////////////////////////////
305                     // The parameter to the Fjt function
306                     const SIMINT_DBLTYPE F_x = SIMINT_MUL(R2, alpha);
307 
308 
309                     const SIMINT_DBLTYPE Q_prefac = mask_load(nlane, Q.prefac + j);
310 
311 
312                     boys_F_split(PRIM_INT__s_s_s_s, F_x, 16);
313                     SIMINT_DBLTYPE prefac = SIMINT_SQRT(one_over_PQalpha_sum);
314                     prefac = SIMINT_MUL(SIMINT_MUL(P_prefac, Q_prefac), prefac);
315                     for(n = 0; n <= 16; n++)
316                         PRIM_INT__s_s_s_s[n] = SIMINT_MUL(PRIM_INT__s_s_s_s[n], prefac);
317 
318                     //////////////////////////////////////////////
319                     // Primitive integrals: Vertical recurrance
320                     //////////////////////////////////////////////
321 
322                     const SIMINT_DBLTYPE vrr_const_1_over_2p = one_over_2p;
323                     const SIMINT_DBLTYPE vrr_const_2_over_2p = SIMINT_MUL(const_2, one_over_2p);
324                     const SIMINT_DBLTYPE vrr_const_3_over_2p = SIMINT_MUL(const_3, one_over_2p);
325                     const SIMINT_DBLTYPE vrr_const_4_over_2p = SIMINT_MUL(const_4, one_over_2p);
326                     const SIMINT_DBLTYPE vrr_const_1_over_2q = one_over_2q;
327                     const SIMINT_DBLTYPE vrr_const_2_over_2q = SIMINT_MUL(const_2, one_over_2q);
328                     const SIMINT_DBLTYPE vrr_const_3_over_2q = SIMINT_MUL(const_3, one_over_2q);
329                     const SIMINT_DBLTYPE vrr_const_4_over_2q = SIMINT_MUL(const_4, one_over_2q);
330                     const SIMINT_DBLTYPE vrr_const_5_over_2q = SIMINT_MUL(const_5, one_over_2q);
331                     const SIMINT_DBLTYPE vrr_const_6_over_2q = SIMINT_MUL(const_6, one_over_2q);
332                     const SIMINT_DBLTYPE vrr_const_7_over_2q = SIMINT_MUL(const_7, one_over_2q);
333                     const SIMINT_DBLTYPE vrr_const_8_over_2q = SIMINT_MUL(const_8, one_over_2q);
334                     const SIMINT_DBLTYPE vrr_const_9_over_2q = SIMINT_MUL(const_9, one_over_2q);
335                     const SIMINT_DBLTYPE vrr_const_10_over_2q = SIMINT_MUL(const_10, one_over_2q);
336                     const SIMINT_DBLTYPE vrr_const_1_over_2pq = one_over_2pq;
337                     const SIMINT_DBLTYPE vrr_const_2_over_2pq = SIMINT_MUL(const_2, one_over_2pq);
338                     const SIMINT_DBLTYPE vrr_const_3_over_2pq = SIMINT_MUL(const_3, one_over_2pq);
339                     const SIMINT_DBLTYPE vrr_const_4_over_2pq = SIMINT_MUL(const_4, one_over_2pq);
340                     const SIMINT_DBLTYPE vrr_const_5_over_2pq = SIMINT_MUL(const_5, one_over_2pq);
341                     const SIMINT_DBLTYPE vrr_const_6_over_2pq = SIMINT_MUL(const_6, one_over_2pq);
342                     const SIMINT_DBLTYPE vrr_const_7_over_2pq = SIMINT_MUL(const_7, one_over_2pq);
343                     const SIMINT_DBLTYPE vrr_const_8_over_2pq = SIMINT_MUL(const_8, one_over_2pq);
344                     const SIMINT_DBLTYPE vrr_const_9_over_2pq = SIMINT_MUL(const_9, one_over_2pq);
345                     const SIMINT_DBLTYPE vrr_const_10_over_2pq = SIMINT_MUL(const_10, one_over_2pq);
346                     const SIMINT_DBLTYPE vrr_const_11_over_2pq = SIMINT_MUL(const_11, one_over_2pq);
347 
348 
349 
350                     // Forming PRIM_INT__s_s_p_s[16 * 3];
351                     for(n = 0; n < 16; ++n)  // loop over orders of auxiliary function
352                     {
353 
354                         PRIM_INT__s_s_p_s[n * 3 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
355                         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]);
356 
357                         PRIM_INT__s_s_p_s[n * 3 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
358                         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]);
359 
360                         PRIM_INT__s_s_p_s[n * 3 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
361                         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]);
362 
363                     }
364 
365 
366 
367                     // Forming PRIM_INT__s_s_d_s[15 * 6];
368                     for(n = 0; n < 15; ++n)  // loop over orders of auxiliary function
369                     {
370 
371                         PRIM_INT__s_s_d_s[n * 6 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_p_s[n * 3 + 0]);
372                         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]);
373                         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]);
374 
375                         PRIM_INT__s_s_d_s[n * 6 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_p_s[n * 3 + 0]);
376                         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]);
377 
378                         PRIM_INT__s_s_d_s[n * 6 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 0]);
379                         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]);
380 
381                         PRIM_INT__s_s_d_s[n * 6 + 3] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_p_s[n * 3 + 1]);
382                         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]);
383                         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]);
384 
385                         PRIM_INT__s_s_d_s[n * 6 + 4] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 1]);
386                         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]);
387 
388                         PRIM_INT__s_s_d_s[n * 6 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_p_s[n * 3 + 2]);
389                         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]);
390                         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]);
391 
392                     }
393 
394 
395 
396                     // Forming PRIM_INT__s_s_f_s[14 * 10];
397                     for(n = 0; n < 14; ++n)  // loop over orders of auxiliary function
398                     {
399 
400                         PRIM_INT__s_s_f_s[n * 10 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 0]);
401                         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]);
402                         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]);
403 
404                         PRIM_INT__s_s_f_s[n * 10 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 0]);
405                         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]);
406 
407                         PRIM_INT__s_s_f_s[n * 10 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 0]);
408                         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]);
409 
410                         PRIM_INT__s_s_f_s[n * 10 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 3]);
411                         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]);
412 
413                         PRIM_INT__s_s_f_s[n * 10 + 4] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 1]);
414                         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]);
415 
416                         PRIM_INT__s_s_f_s[n * 10 + 5] = SIMINT_MUL(Q_PA[0], PRIM_INT__s_s_d_s[n * 6 + 5]);
417                         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]);
418 
419                         PRIM_INT__s_s_f_s[n * 10 + 6] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 3]);
420                         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]);
421                         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]);
422 
423                         PRIM_INT__s_s_f_s[n * 10 + 7] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 3]);
424                         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]);
425 
426                         PRIM_INT__s_s_f_s[n * 10 + 8] = SIMINT_MUL(Q_PA[1], PRIM_INT__s_s_d_s[n * 6 + 5]);
427                         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]);
428 
429                         PRIM_INT__s_s_f_s[n * 10 + 9] = SIMINT_MUL(Q_PA[2], PRIM_INT__s_s_d_s[n * 6 + 5]);
430                         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]);
431                         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]);
432 
433                     }
434 
435 
436                     VRR_K_s_s_g_s(
437                             PRIM_INT__s_s_g_s,
438                             PRIM_INT__s_s_f_s,
439                             PRIM_INT__s_s_d_s,
440                             Q_PA,
441                             a_over_q,
442                             aoq_PQ,
443                             one_over_2q,
444                             13);
445 
446 
447                     VRR_K_s_s_h_s(
448                             PRIM_INT__s_s_h_s,
449                             PRIM_INT__s_s_g_s,
450                             PRIM_INT__s_s_f_s,
451                             Q_PA,
452                             a_over_q,
453                             aoq_PQ,
454                             one_over_2q,
455                             12);
456 
457 
458                     ostei_general_vrr1_K(6, 11,
459                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
460                             PRIM_INT__s_s_h_s, PRIM_INT__s_s_g_s, PRIM_INT__s_s_i_s);
461 
462 
463                     ostei_general_vrr_I(1, 0, 6, 0, 5,
464                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
465                             PRIM_INT__s_s_i_s, NULL, NULL, PRIM_INT__s_s_h_s, NULL, PRIM_INT__p_s_i_s);
466 
467 
468                     ostei_general_vrr_I(1, 0, 5, 0, 5,
469                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
470                             PRIM_INT__s_s_h_s, NULL, NULL, PRIM_INT__s_s_g_s, NULL, PRIM_INT__p_s_h_s);
471 
472 
473                     ostei_general_vrr_I(2, 0, 6, 0, 4,
474                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
475                             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);
476 
477 
478                     VRR_I_p_s_g_s(
479                             PRIM_INT__p_s_g_s,
480                             PRIM_INT__s_s_g_s,
481                             PRIM_INT__s_s_f_s,
482                             P_PA,
483                             aop_PQ,
484                             one_over_2pq,
485                             5);
486 
487 
488                     ostei_general_vrr_I(2, 0, 5, 0, 4,
489                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
490                             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);
491 
492 
493                     ostei_general_vrr_I(3, 0, 6, 0, 3,
494                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
495                             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);
496 
497 
498                     VRR_I_p_s_f_s(
499                             PRIM_INT__p_s_f_s,
500                             PRIM_INT__s_s_f_s,
501                             PRIM_INT__s_s_d_s,
502                             P_PA,
503                             aop_PQ,
504                             one_over_2pq,
505                             5);
506 
507 
508                     ostei_general_vrr_I(2, 0, 4, 0, 4,
509                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
510                             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);
511 
512 
513                     ostei_general_vrr_I(3, 0, 5, 0, 3,
514                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
515                             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);
516 
517 
518                     ostei_general_vrr_I(4, 0, 6, 0, 2,
519                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
520                             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);
521 
522 
523                     ostei_general_vrr1_K(7, 10,
524                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
525                             PRIM_INT__s_s_i_s, PRIM_INT__s_s_h_s, PRIM_INT__s_s_k_s);
526 
527 
528                     ostei_general_vrr_I(1, 0, 7, 0, 5,
529                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
530                             PRIM_INT__s_s_k_s, NULL, NULL, PRIM_INT__s_s_i_s, NULL, PRIM_INT__p_s_k_s);
531 
532 
533                     ostei_general_vrr_I(2, 0, 7, 0, 4,
534                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
535                             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);
536 
537 
538                     ostei_general_vrr_I(3, 0, 7, 0, 3,
539                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
540                             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);
541 
542 
543                     ostei_general_vrr_I(4, 0, 7, 0, 2,
544                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
545                             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);
546 
547 
548 
549                     // Forming PRIM_INT__p_s_d_s[5 * 18];
550                     for(n = 0; n < 5; ++n)  // loop over orders of auxiliary function
551                     {
552 
553                         PRIM_INT__p_s_d_s[n * 18 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_d_s[n * 6 + 0]);
554                         PRIM_INT__p_s_d_s[n * 18 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__p_s_d_s[n * 18 + 0]);
555                         PRIM_INT__p_s_d_s[n * 18 + 0] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__s_s_p_s[(n+1) * 3 + 0], PRIM_INT__p_s_d_s[n * 18 + 0]);
556 
557                         PRIM_INT__p_s_d_s[n * 18 + 1] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_d_s[n * 6 + 1]);
558                         PRIM_INT__p_s_d_s[n * 18 + 1] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 1], PRIM_INT__p_s_d_s[n * 18 + 1]);
559                         PRIM_INT__p_s_d_s[n * 18 + 1] = 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 + 1]);
560 
561                         PRIM_INT__p_s_d_s[n * 18 + 2] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_d_s[n * 6 + 2]);
562                         PRIM_INT__p_s_d_s[n * 18 + 2] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 2], PRIM_INT__p_s_d_s[n * 18 + 2]);
563                         PRIM_INT__p_s_d_s[n * 18 + 2] = 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 + 2]);
564 
565                         PRIM_INT__p_s_d_s[n * 18 + 3] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_d_s[n * 6 + 3]);
566                         PRIM_INT__p_s_d_s[n * 18 + 3] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__p_s_d_s[n * 18 + 3]);
567 
568                         PRIM_INT__p_s_d_s[n * 18 + 4] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_d_s[n * 6 + 4]);
569                         PRIM_INT__p_s_d_s[n * 18 + 4] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 4], PRIM_INT__p_s_d_s[n * 18 + 4]);
570 
571                         PRIM_INT__p_s_d_s[n * 18 + 5] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_d_s[n * 6 + 5]);
572                         PRIM_INT__p_s_d_s[n * 18 + 5] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__p_s_d_s[n * 18 + 5]);
573 
574                         PRIM_INT__p_s_d_s[n * 18 + 6] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_d_s[n * 6 + 0]);
575                         PRIM_INT__p_s_d_s[n * 18 + 6] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__p_s_d_s[n * 18 + 6]);
576 
577                         PRIM_INT__p_s_d_s[n * 18 + 7] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_d_s[n * 6 + 1]);
578                         PRIM_INT__p_s_d_s[n * 18 + 7] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 1], PRIM_INT__p_s_d_s[n * 18 + 7]);
579                         PRIM_INT__p_s_d_s[n * 18 + 7] = 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 + 7]);
580 
581                         PRIM_INT__p_s_d_s[n * 18 + 8] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_d_s[n * 6 + 2]);
582                         PRIM_INT__p_s_d_s[n * 18 + 8] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 2], PRIM_INT__p_s_d_s[n * 18 + 8]);
583 
584                         PRIM_INT__p_s_d_s[n * 18 + 9] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_d_s[n * 6 + 3]);
585                         PRIM_INT__p_s_d_s[n * 18 + 9] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__p_s_d_s[n * 18 + 9]);
586                         PRIM_INT__p_s_d_s[n * 18 + 9] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__s_s_p_s[(n+1) * 3 + 1], PRIM_INT__p_s_d_s[n * 18 + 9]);
587 
588                         PRIM_INT__p_s_d_s[n * 18 + 10] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_d_s[n * 6 + 4]);
589                         PRIM_INT__p_s_d_s[n * 18 + 10] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 4], PRIM_INT__p_s_d_s[n * 18 + 10]);
590                         PRIM_INT__p_s_d_s[n * 18 + 10] = 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 + 10]);
591 
592                         PRIM_INT__p_s_d_s[n * 18 + 11] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_d_s[n * 6 + 5]);
593                         PRIM_INT__p_s_d_s[n * 18 + 11] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__p_s_d_s[n * 18 + 11]);
594 
595                         PRIM_INT__p_s_d_s[n * 18 + 12] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_d_s[n * 6 + 0]);
596                         PRIM_INT__p_s_d_s[n * 18 + 12] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 0], PRIM_INT__p_s_d_s[n * 18 + 12]);
597 
598                         PRIM_INT__p_s_d_s[n * 18 + 13] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_d_s[n * 6 + 1]);
599                         PRIM_INT__p_s_d_s[n * 18 + 13] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 1], PRIM_INT__p_s_d_s[n * 18 + 13]);
600 
601                         PRIM_INT__p_s_d_s[n * 18 + 14] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_d_s[n * 6 + 2]);
602                         PRIM_INT__p_s_d_s[n * 18 + 14] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 2], PRIM_INT__p_s_d_s[n * 18 + 14]);
603                         PRIM_INT__p_s_d_s[n * 18 + 14] = 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 + 14]);
604 
605                         PRIM_INT__p_s_d_s[n * 18 + 15] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_d_s[n * 6 + 3]);
606                         PRIM_INT__p_s_d_s[n * 18 + 15] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 3], PRIM_INT__p_s_d_s[n * 18 + 15]);
607 
608                         PRIM_INT__p_s_d_s[n * 18 + 16] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_d_s[n * 6 + 4]);
609                         PRIM_INT__p_s_d_s[n * 18 + 16] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 4], PRIM_INT__p_s_d_s[n * 18 + 16]);
610                         PRIM_INT__p_s_d_s[n * 18 + 16] = 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 + 16]);
611 
612                         PRIM_INT__p_s_d_s[n * 18 + 17] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_d_s[n * 6 + 5]);
613                         PRIM_INT__p_s_d_s[n * 18 + 17] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_d_s[(n+1) * 6 + 5], PRIM_INT__p_s_d_s[n * 18 + 17]);
614                         PRIM_INT__p_s_d_s[n * 18 + 17] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__s_s_p_s[(n+1) * 3 + 2], PRIM_INT__p_s_d_s[n * 18 + 17]);
615 
616                     }
617 
618 
619                     VRR_I_d_s_f_s(
620                             PRIM_INT__d_s_f_s,
621                             PRIM_INT__p_s_f_s,
622                             PRIM_INT__s_s_f_s,
623                             PRIM_INT__p_s_d_s,
624                             P_PA,
625                             a_over_p,
626                             aop_PQ,
627                             one_over_2p,
628                             one_over_2pq,
629                             4);
630 
631 
632                     ostei_general_vrr_I(3, 0, 4, 0, 3,
633                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
634                             PRIM_INT__d_s_g_s, PRIM_INT__p_s_g_s, NULL, PRIM_INT__d_s_f_s, NULL, PRIM_INT__f_s_g_s);
635 
636 
637                     ostei_general_vrr_I(4, 0, 5, 0, 2,
638                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
639                             PRIM_INT__f_s_h_s, PRIM_INT__d_s_h_s, NULL, PRIM_INT__f_s_g_s, NULL, PRIM_INT__g_s_h_s);
640 
641 
642                     ostei_general_vrr_I(5, 0, 6, 0, 1,
643                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
644                             PRIM_INT__g_s_i_s, PRIM_INT__f_s_i_s, NULL, PRIM_INT__g_s_h_s, NULL, PRIM_INT__h_s_i_s);
645 
646 
647                     ostei_general_vrr1_K(8, 9,
648                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
649                             PRIM_INT__s_s_k_s, PRIM_INT__s_s_i_s, PRIM_INT__s_s_l_s);
650 
651 
652                     ostei_general_vrr_I(1, 0, 8, 0, 5,
653                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
654                             PRIM_INT__s_s_l_s, NULL, NULL, PRIM_INT__s_s_k_s, NULL, PRIM_INT__p_s_l_s);
655 
656 
657                     ostei_general_vrr_I(2, 0, 8, 0, 4,
658                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
659                             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);
660 
661 
662                     ostei_general_vrr_I(3, 0, 8, 0, 3,
663                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
664                             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);
665 
666 
667                     ostei_general_vrr_I(4, 0, 8, 0, 2,
668                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
669                             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);
670 
671 
672                     ostei_general_vrr_I(5, 0, 7, 0, 1,
673                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
674                             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);
675 
676 
677                     ostei_general_vrr1_K(9, 8,
678                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
679                             PRIM_INT__s_s_l_s, PRIM_INT__s_s_k_s, PRIM_INT__s_s_m_s);
680 
681 
682                     ostei_general_vrr_I(1, 0, 9, 0, 5,
683                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
684                             PRIM_INT__s_s_m_s, NULL, NULL, PRIM_INT__s_s_l_s, NULL, PRIM_INT__p_s_m_s);
685 
686 
687                     ostei_general_vrr_I(2, 0, 9, 0, 4,
688                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
689                             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);
690 
691 
692                     ostei_general_vrr_I(3, 0, 9, 0, 3,
693                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
694                             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);
695 
696 
697                     ostei_general_vrr_I(4, 0, 9, 0, 2,
698                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
699                             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);
700 
701 
702                     ostei_general_vrr_I(5, 0, 8, 0, 1,
703                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
704                             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);
705 
706 
707                     ostei_general_vrr1_K(10, 7,
708                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
709                             PRIM_INT__s_s_m_s, PRIM_INT__s_s_l_s, PRIM_INT__s_s_n_s);
710 
711 
712                     ostei_general_vrr_I(1, 0, 10, 0, 5,
713                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
714                             PRIM_INT__s_s_n_s, NULL, NULL, PRIM_INT__s_s_m_s, NULL, PRIM_INT__p_s_n_s);
715 
716 
717                     ostei_general_vrr_I(2, 0, 10, 0, 4,
718                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
719                             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);
720 
721 
722                     ostei_general_vrr_I(3, 0, 10, 0, 3,
723                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
724                             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);
725 
726 
727                     ostei_general_vrr_I(4, 0, 10, 0, 2,
728                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
729                             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);
730 
731 
732                     ostei_general_vrr_I(5, 0, 9, 0, 1,
733                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
734                             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);
735 
736 
737                     ostei_general_vrr1_K(11, 6,
738                             one_over_2q, a_over_q, aoq_PQ, Q_PA,
739                             PRIM_INT__s_s_n_s, PRIM_INT__s_s_m_s, PRIM_INT__s_s_o_s);
740 
741 
742                     ostei_general_vrr_I(1, 0, 11, 0, 5,
743                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
744                             PRIM_INT__s_s_o_s, NULL, NULL, PRIM_INT__s_s_n_s, NULL, PRIM_INT__p_s_o_s);
745 
746 
747                     ostei_general_vrr_I(2, 0, 11, 0, 4,
748                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
749                             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);
750 
751 
752                     ostei_general_vrr_I(3, 0, 11, 0, 3,
753                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
754                             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);
755 
756 
757                     ostei_general_vrr_I(4, 0, 11, 0, 2,
758                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
759                             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);
760 
761 
762                     ostei_general_vrr_I(5, 0, 10, 0, 1,
763                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
764                             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);
765 
766 
767                     ostei_general_vrr_I(5, 0, 11, 0, 1,
768                             one_over_2p, a_over_p, one_over_2pq, aop_PQ, P_PA,
769                             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);
770 
771 
772 
773 
774                     ////////////////////////////////////
775                     // Accumulate contracted integrals
776                     ////////////////////////////////////
777                     if(lastoffset == 0)
778                     {
779                         contract_all(420, PRIM_INT__g_s_i_s, PRIM_PTR_INT__g_s_i_s);
780                         contract_all(540, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
781                         contract_all(675, PRIM_INT__g_s_l_s, PRIM_PTR_INT__g_s_l_s);
782                         contract_all(825, PRIM_INT__g_s_m_s, PRIM_PTR_INT__g_s_m_s);
783                         contract_all(990, PRIM_INT__g_s_n_s, PRIM_PTR_INT__g_s_n_s);
784                         contract_all(1170, PRIM_INT__g_s_o_s, PRIM_PTR_INT__g_s_o_s);
785                         contract_all(588, PRIM_INT__h_s_i_s, PRIM_PTR_INT__h_s_i_s);
786                         contract_all(756, PRIM_INT__h_s_k_s, PRIM_PTR_INT__h_s_k_s);
787                         contract_all(945, PRIM_INT__h_s_l_s, PRIM_PTR_INT__h_s_l_s);
788                         contract_all(1155, PRIM_INT__h_s_m_s, PRIM_PTR_INT__h_s_m_s);
789                         contract_all(1386, PRIM_INT__h_s_n_s, PRIM_PTR_INT__h_s_n_s);
790                         contract_all(1638, PRIM_INT__h_s_o_s, PRIM_PTR_INT__h_s_o_s);
791                     }
792                     else
793                     {
794                         contract(420, shelloffsets, PRIM_INT__g_s_i_s, PRIM_PTR_INT__g_s_i_s);
795                         contract(540, shelloffsets, PRIM_INT__g_s_k_s, PRIM_PTR_INT__g_s_k_s);
796                         contract(675, shelloffsets, PRIM_INT__g_s_l_s, PRIM_PTR_INT__g_s_l_s);
797                         contract(825, shelloffsets, PRIM_INT__g_s_m_s, PRIM_PTR_INT__g_s_m_s);
798                         contract(990, shelloffsets, PRIM_INT__g_s_n_s, PRIM_PTR_INT__g_s_n_s);
799                         contract(1170, shelloffsets, PRIM_INT__g_s_o_s, PRIM_PTR_INT__g_s_o_s);
800                         contract(588, shelloffsets, PRIM_INT__h_s_i_s, PRIM_PTR_INT__h_s_i_s);
801                         contract(756, shelloffsets, PRIM_INT__h_s_k_s, PRIM_PTR_INT__h_s_k_s);
802                         contract(945, shelloffsets, PRIM_INT__h_s_l_s, PRIM_PTR_INT__h_s_l_s);
803                         contract(1155, shelloffsets, PRIM_INT__h_s_m_s, PRIM_PTR_INT__h_s_m_s);
804                         contract(1386, shelloffsets, PRIM_INT__h_s_n_s, PRIM_PTR_INT__h_s_n_s);
805                         contract(1638, shelloffsets, PRIM_INT__h_s_o_s, PRIM_PTR_INT__h_s_o_s);
806                         PRIM_PTR_INT__g_s_i_s += lastoffset*420;
807                         PRIM_PTR_INT__g_s_k_s += lastoffset*540;
808                         PRIM_PTR_INT__g_s_l_s += lastoffset*675;
809                         PRIM_PTR_INT__g_s_m_s += lastoffset*825;
810                         PRIM_PTR_INT__g_s_n_s += lastoffset*990;
811                         PRIM_PTR_INT__g_s_o_s += lastoffset*1170;
812                         PRIM_PTR_INT__h_s_i_s += lastoffset*588;
813                         PRIM_PTR_INT__h_s_k_s += lastoffset*756;
814                         PRIM_PTR_INT__h_s_l_s += lastoffset*945;
815                         PRIM_PTR_INT__h_s_m_s += lastoffset*1155;
816                         PRIM_PTR_INT__h_s_n_s += lastoffset*1386;
817                         PRIM_PTR_INT__h_s_o_s += lastoffset*1638;
818                     }
819 
820                 }  // close loop over j
821             }  // close loop over i
822 
823             //Advance to the next batch
824             jstart = SIMINT_SIMD_ROUND(jend);
825 
826             //////////////////////////////////////////////
827             // Contracted integrals: Horizontal recurrance
828             //////////////////////////////////////////////
829 
830 
831             const double hAB[3] = { P.AB_x[ab], P.AB_y[ab], P.AB_z[ab] };
832 
833 
834             for(abcd = 0; abcd < nshellbatch; ++abcd, ++real_abcd)
835             {
836                 const double hCD[3] = { Q.AB_x[cd+abcd], Q.AB_y[cd+abcd], Q.AB_z[cd+abcd] };
837 
838                 // set up HRR pointers
839                 double const * restrict HRR_INT__g_s_i_s = INT__g_s_i_s + abcd * 420;
840                 double const * restrict HRR_INT__g_s_k_s = INT__g_s_k_s + abcd * 540;
841                 double const * restrict HRR_INT__g_s_l_s = INT__g_s_l_s + abcd * 675;
842                 double const * restrict HRR_INT__g_s_m_s = INT__g_s_m_s + abcd * 825;
843                 double const * restrict HRR_INT__g_s_n_s = INT__g_s_n_s + abcd * 990;
844                 double const * restrict HRR_INT__g_s_o_s = INT__g_s_o_s + abcd * 1170;
845                 double const * restrict HRR_INT__h_s_i_s = INT__h_s_i_s + abcd * 588;
846                 double const * restrict HRR_INT__h_s_k_s = INT__h_s_k_s + abcd * 756;
847                 double const * restrict HRR_INT__h_s_l_s = INT__h_s_l_s + abcd * 945;
848                 double const * restrict HRR_INT__h_s_m_s = INT__h_s_m_s + abcd * 1155;
849                 double const * restrict HRR_INT__h_s_n_s = INT__h_s_n_s + abcd * 1386;
850                 double const * restrict HRR_INT__h_s_o_s = INT__h_s_o_s + abcd * 1638;
851                 double * restrict HRR_INT__g_p_i_h = INT__g_p_i_h + real_abcd * 26460;
852 
853                 // form INT__g_p_i_s
854                 HRR_J_g_p(
855                     HRR_INT__g_p_i_s,
856                     HRR_INT__g_s_i_s,
857                     HRR_INT__h_s_i_s,
858                     hAB, 28);
859 
860                 // form INT__g_p_k_s
861                 HRR_J_g_p(
862                     HRR_INT__g_p_k_s,
863                     HRR_INT__g_s_k_s,
864                     HRR_INT__h_s_k_s,
865                     hAB, 36);
866 
867                 // form INT__g_p_l_s
868                 HRR_J_g_p(
869                     HRR_INT__g_p_l_s,
870                     HRR_INT__g_s_l_s,
871                     HRR_INT__h_s_l_s,
872                     hAB, 45);
873 
874                 // form INT__g_p_m_s
875                 HRR_J_g_p(
876                     HRR_INT__g_p_m_s,
877                     HRR_INT__g_s_m_s,
878                     HRR_INT__h_s_m_s,
879                     hAB, 55);
880 
881                 // form INT__g_p_n_s
882                 HRR_J_g_p(
883                     HRR_INT__g_p_n_s,
884                     HRR_INT__g_s_n_s,
885                     HRR_INT__h_s_n_s,
886                     hAB, 66);
887 
888                 // form INT__g_p_o_s
889                 HRR_J_g_p(
890                     HRR_INT__g_p_o_s,
891                     HRR_INT__g_s_o_s,
892                     HRR_INT__h_s_o_s,
893                     hAB, 78);
894 
895                 // form INT__g_p_i_p
896                 ostei_general_hrr_L(4, 1, 6, 1, hCD, HRR_INT__g_p_k_s, HRR_INT__g_p_i_s, HRR_INT__g_p_i_p);
897 
898                 // form INT__g_p_k_p
899                 ostei_general_hrr_L(4, 1, 7, 1, hCD, HRR_INT__g_p_l_s, HRR_INT__g_p_k_s, HRR_INT__g_p_k_p);
900 
901                 // form INT__g_p_l_p
902                 ostei_general_hrr_L(4, 1, 8, 1, hCD, HRR_INT__g_p_m_s, HRR_INT__g_p_l_s, HRR_INT__g_p_l_p);
903 
904                 // form INT__g_p_m_p
905                 ostei_general_hrr_L(4, 1, 9, 1, hCD, HRR_INT__g_p_n_s, HRR_INT__g_p_m_s, HRR_INT__g_p_m_p);
906 
907                 // form INT__g_p_n_p
908                 ostei_general_hrr_L(4, 1, 10, 1, hCD, HRR_INT__g_p_o_s, HRR_INT__g_p_n_s, HRR_INT__g_p_n_p);
909 
910                 // form INT__g_p_i_d
911                 ostei_general_hrr_L(4, 1, 6, 2, hCD, HRR_INT__g_p_k_p, HRR_INT__g_p_i_p, HRR_INT__g_p_i_d);
912 
913                 // form INT__g_p_k_d
914                 ostei_general_hrr_L(4, 1, 7, 2, hCD, HRR_INT__g_p_l_p, HRR_INT__g_p_k_p, HRR_INT__g_p_k_d);
915 
916                 // form INT__g_p_l_d
917                 ostei_general_hrr_L(4, 1, 8, 2, hCD, HRR_INT__g_p_m_p, HRR_INT__g_p_l_p, HRR_INT__g_p_l_d);
918 
919                 // form INT__g_p_m_d
920                 ostei_general_hrr_L(4, 1, 9, 2, hCD, HRR_INT__g_p_n_p, HRR_INT__g_p_m_p, HRR_INT__g_p_m_d);
921 
922                 // form INT__g_p_i_f
923                 ostei_general_hrr_L(4, 1, 6, 3, hCD, HRR_INT__g_p_k_d, HRR_INT__g_p_i_d, HRR_INT__g_p_i_f);
924 
925                 // form INT__g_p_k_f
926                 ostei_general_hrr_L(4, 1, 7, 3, hCD, HRR_INT__g_p_l_d, HRR_INT__g_p_k_d, HRR_INT__g_p_k_f);
927 
928                 // form INT__g_p_l_f
929                 ostei_general_hrr_L(4, 1, 8, 3, hCD, HRR_INT__g_p_m_d, HRR_INT__g_p_l_d, HRR_INT__g_p_l_f);
930 
931                 // form INT__g_p_i_g
932                 ostei_general_hrr_L(4, 1, 6, 4, hCD, HRR_INT__g_p_k_f, HRR_INT__g_p_i_f, HRR_INT__g_p_i_g);
933 
934                 // form INT__g_p_k_g
935                 ostei_general_hrr_L(4, 1, 7, 4, hCD, HRR_INT__g_p_l_f, HRR_INT__g_p_k_f, HRR_INT__g_p_k_g);
936 
937                 // form INT__g_p_i_h
938                 ostei_general_hrr_L(4, 1, 6, 5, hCD, HRR_INT__g_p_k_g, HRR_INT__g_p_i_g, HRR_INT__g_p_i_h);
939 
940 
941             }  // close HRR loop
942 
943 
944         }   // close loop cdbatch
945 
946         istart = iend;
947     }  // close loop over ab
948 
949     return P.nshell12_clip * Q.nshell12_clip;
950 }
951 
ostei_p_g_i_h(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__p_g_i_h)952 int ostei_p_g_i_h(struct simint_multi_shellpair const P,
953                   struct simint_multi_shellpair const Q,
954                   double screen_tol,
955                   double * const restrict work,
956                   double * const restrict INT__p_g_i_h)
957 {
958     double P_AB[3*P.nshell12];
959     struct simint_multi_shellpair P_tmp = P;
960     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
961     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
962     P_tmp.AB_x = P_AB;
963     P_tmp.AB_y = P_AB + P.nshell12;
964     P_tmp.AB_z = P_AB + 2*P.nshell12;
965 
966     for(int i = 0; i < P.nshell12; i++)
967     {
968         P_tmp.AB_x[i] = -P.AB_x[i];
969         P_tmp.AB_y[i] = -P.AB_y[i];
970         P_tmp.AB_z[i] = -P.AB_z[i];
971     }
972 
973     int ret = ostei_g_p_i_h(P_tmp, Q, screen_tol, work, INT__p_g_i_h);
974     double buffer[26460] SIMINT_ALIGN_ARRAY_DBL;
975 
976     for(int q = 0; q < ret; q++)
977     {
978         int idx = 0;
979         for(int a = 0; a < 3; ++a)
980         for(int b = 0; b < 15; ++b)
981         for(int c = 0; c < 28; ++c)
982         for(int d = 0; d < 21; ++d)
983             buffer[idx++] = INT__p_g_i_h[q*26460+b*1764+a*588+c*21+d];
984 
985         memcpy(INT__p_g_i_h+q*26460, buffer, 26460*sizeof(double));
986     }
987 
988     return ret;
989 }
990 
ostei_g_p_h_i(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__g_p_h_i)991 int ostei_g_p_h_i(struct simint_multi_shellpair const P,
992                   struct simint_multi_shellpair const Q,
993                   double screen_tol,
994                   double * const restrict work,
995                   double * const restrict INT__g_p_h_i)
996 {
997     double Q_AB[3*Q.nshell12];
998     struct simint_multi_shellpair Q_tmp = Q;
999     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
1000     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
1001     Q_tmp.AB_x = Q_AB;
1002     Q_tmp.AB_y = Q_AB + Q.nshell12;
1003     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
1004 
1005     for(int i = 0; i < Q.nshell12; i++)
1006     {
1007         Q_tmp.AB_x[i] = -Q.AB_x[i];
1008         Q_tmp.AB_y[i] = -Q.AB_y[i];
1009         Q_tmp.AB_z[i] = -Q.AB_z[i];
1010     }
1011 
1012     int ret = ostei_g_p_i_h(P, Q_tmp, screen_tol, work, INT__g_p_h_i);
1013     double buffer[26460] SIMINT_ALIGN_ARRAY_DBL;
1014 
1015     for(int q = 0; q < ret; q++)
1016     {
1017         int idx = 0;
1018         for(int a = 0; a < 15; ++a)
1019         for(int b = 0; b < 3; ++b)
1020         for(int c = 0; c < 21; ++c)
1021         for(int d = 0; d < 28; ++d)
1022             buffer[idx++] = INT__g_p_h_i[q*26460+a*1764+b*588+d*21+c];
1023 
1024         memcpy(INT__g_p_h_i+q*26460, buffer, 26460*sizeof(double));
1025     }
1026 
1027     return ret;
1028 }
1029 
ostei_p_g_h_i(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__p_g_h_i)1030 int ostei_p_g_h_i(struct simint_multi_shellpair const P,
1031                   struct simint_multi_shellpair const Q,
1032                   double screen_tol,
1033                   double * const restrict work,
1034                   double * const restrict INT__p_g_h_i)
1035 {
1036     double P_AB[3*P.nshell12];
1037     struct simint_multi_shellpair P_tmp = P;
1038     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
1039     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
1040     P_tmp.AB_x = P_AB;
1041     P_tmp.AB_y = P_AB + P.nshell12;
1042     P_tmp.AB_z = P_AB + 2*P.nshell12;
1043 
1044     for(int i = 0; i < P.nshell12; i++)
1045     {
1046         P_tmp.AB_x[i] = -P.AB_x[i];
1047         P_tmp.AB_y[i] = -P.AB_y[i];
1048         P_tmp.AB_z[i] = -P.AB_z[i];
1049     }
1050 
1051     double Q_AB[3*Q.nshell12];
1052     struct simint_multi_shellpair Q_tmp = Q;
1053     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
1054     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
1055     Q_tmp.AB_x = Q_AB;
1056     Q_tmp.AB_y = Q_AB + Q.nshell12;
1057     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
1058 
1059     for(int i = 0; i < Q.nshell12; i++)
1060     {
1061         Q_tmp.AB_x[i] = -Q.AB_x[i];
1062         Q_tmp.AB_y[i] = -Q.AB_y[i];
1063         Q_tmp.AB_z[i] = -Q.AB_z[i];
1064     }
1065 
1066     int ret = ostei_g_p_i_h(P_tmp, Q_tmp, screen_tol, work, INT__p_g_h_i);
1067     double buffer[26460] SIMINT_ALIGN_ARRAY_DBL;
1068 
1069     for(int q = 0; q < ret; q++)
1070     {
1071         int idx = 0;
1072         for(int a = 0; a < 3; ++a)
1073         for(int b = 0; b < 15; ++b)
1074         for(int c = 0; c < 21; ++c)
1075         for(int d = 0; d < 28; ++d)
1076             buffer[idx++] = INT__p_g_h_i[q*26460+b*1764+a*588+d*21+c];
1077 
1078         memcpy(INT__p_g_h_i+q*26460, buffer, 26460*sizeof(double));
1079     }
1080 
1081     return ret;
1082 }
1083 
1084