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_i_g_g_p(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__i_g_g_p)8 int ostei_i_g_g_p(struct simint_multi_shellpair const P,
9                   struct simint_multi_shellpair const Q,
10                   double screen_tol,
11                   double * const restrict work,
12                   double * const restrict INT__i_g_g_p)
13 {
14 
15     SIMINT_ASSUME_ALIGN_DBL(work);
16     SIMINT_ASSUME_ALIGN_DBL(INT__i_g_g_p);
17     int ab, cd, abcd;
18     int istart, jstart;
19     int iprimcd, nprim_icd, icd;
20     const int check_screen = (screen_tol > 0.0);
21     int i, j;
22     int n;
23     int not_screened;
24     int real_abcd;
25     int iket;
26     int ibra;
27 
28     // partition workspace
29     double * const INT__i_s_g_s = work + (SIMINT_NSHELL_SIMD * 0);
30     double * const INT__i_s_h_s = work + (SIMINT_NSHELL_SIMD * 420);
31     double * const INT__k_s_g_s = work + (SIMINT_NSHELL_SIMD * 1008);
32     double * const INT__k_s_h_s = work + (SIMINT_NSHELL_SIMD * 1548);
33     double * const INT__l_s_g_s = work + (SIMINT_NSHELL_SIMD * 2304);
34     double * const INT__l_s_h_s = work + (SIMINT_NSHELL_SIMD * 2979);
35     double * const INT__m_s_g_s = work + (SIMINT_NSHELL_SIMD * 3924);
36     double * const INT__m_s_h_s = work + (SIMINT_NSHELL_SIMD * 4749);
37     double * const INT__n_s_g_s = work + (SIMINT_NSHELL_SIMD * 5904);
38     double * const INT__n_s_h_s = work + (SIMINT_NSHELL_SIMD * 6894);
39     SIMINT_DBLTYPE * const primwork = (SIMINT_DBLTYPE *)(work + SIMINT_NSHELL_SIMD*8280);
40     SIMINT_DBLTYPE * const restrict PRIM_INT__s_s_s_s = primwork + 0;
41     SIMINT_DBLTYPE * const restrict PRIM_INT__p_s_s_s = primwork + 16;
42     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_s_s = primwork + 61;
43     SIMINT_DBLTYPE * const restrict PRIM_INT__d_s_p_s = primwork + 145;
44     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_s_s = primwork + 235;
45     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_p_s = primwork + 365;
46     SIMINT_DBLTYPE * const restrict PRIM_INT__f_s_d_s = primwork + 515;
47     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_s_s = primwork + 755;
48     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_p_s = primwork + 935;
49     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_d_s = primwork + 1160;
50     SIMINT_DBLTYPE * const restrict PRIM_INT__g_s_f_s = primwork + 1520;
51     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_s_s = primwork + 1970;
52     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_p_s = primwork + 2201;
53     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_d_s = primwork + 2516;
54     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_f_s = primwork + 3020;
55     SIMINT_DBLTYPE * const restrict PRIM_INT__h_s_g_s = primwork + 3650;
56     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_s_s = primwork + 4280;
57     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_p_s = primwork + 4560;
58     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_d_s = primwork + 4980;
59     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_f_s = primwork + 5652;
60     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_g_s = primwork + 6492;
61     SIMINT_DBLTYPE * const restrict PRIM_INT__i_s_h_s = primwork + 7332;
62     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_s_s = primwork + 7920;
63     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_p_s = primwork + 8244;
64     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_d_s = primwork + 8784;
65     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_f_s = primwork + 9648;
66     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_g_s = primwork + 10728;
67     SIMINT_DBLTYPE * const restrict PRIM_INT__k_s_h_s = primwork + 11808;
68     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_s_s = primwork + 12564;
69     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_p_s = primwork + 12924;
70     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_d_s = primwork + 13599;
71     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_f_s = primwork + 14679;
72     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_g_s = primwork + 16029;
73     SIMINT_DBLTYPE * const restrict PRIM_INT__l_s_h_s = primwork + 17379;
74     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_s_s = primwork + 18324;
75     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_p_s = primwork + 18709;
76     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_d_s = primwork + 19534;
77     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_f_s = primwork + 20854;
78     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_g_s = primwork + 22504;
79     SIMINT_DBLTYPE * const restrict PRIM_INT__m_s_h_s = primwork + 24154;
80     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_s_s = primwork + 25309;
81     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_p_s = primwork + 25705;
82     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_d_s = primwork + 26695;
83     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_f_s = primwork + 28279;
84     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_g_s = primwork + 30259;
85     SIMINT_DBLTYPE * const restrict PRIM_INT__n_s_h_s = primwork + 32239;
86     double * const hrrwork = (double *)(primwork + 33625);
87     double * const HRR_INT__i_p_g_s = hrrwork + 0;
88     double * const HRR_INT__i_p_h_s = hrrwork + 1260;
89     double * const HRR_INT__i_d_g_s = hrrwork + 3024;
90     double * const HRR_INT__i_d_h_s = hrrwork + 5544;
91     double * const HRR_INT__i_f_g_s = hrrwork + 9072;
92     double * const HRR_INT__i_f_h_s = hrrwork + 13272;
93     double * const HRR_INT__i_g_g_s = hrrwork + 19152;
94     double * const HRR_INT__i_g_h_s = hrrwork + 25452;
95     double * const HRR_INT__k_p_g_s = hrrwork + 34272;
96     double * const HRR_INT__k_p_h_s = hrrwork + 35892;
97     double * const HRR_INT__k_d_g_s = hrrwork + 38160;
98     double * const HRR_INT__k_d_h_s = hrrwork + 41400;
99     double * const HRR_INT__k_f_g_s = hrrwork + 45936;
100     double * const HRR_INT__k_f_h_s = hrrwork + 51336;
101     double * const HRR_INT__l_p_g_s = hrrwork + 58896;
102     double * const HRR_INT__l_p_h_s = hrrwork + 60921;
103     double * const HRR_INT__l_d_g_s = hrrwork + 63756;
104     double * const HRR_INT__l_d_h_s = hrrwork + 67806;
105     double * const HRR_INT__m_p_g_s = hrrwork + 73476;
106     double * const HRR_INT__m_p_h_s = hrrwork + 75951;
107 
108 
109     // Create constants
110     const SIMINT_DBLTYPE const_1 = SIMINT_DBLSET1(1);
111     const SIMINT_DBLTYPE const_10 = SIMINT_DBLSET1(10);
112     const SIMINT_DBLTYPE const_2 = SIMINT_DBLSET1(2);
113     const SIMINT_DBLTYPE const_3 = SIMINT_DBLSET1(3);
114     const SIMINT_DBLTYPE const_4 = SIMINT_DBLSET1(4);
115     const SIMINT_DBLTYPE const_5 = SIMINT_DBLSET1(5);
116     const SIMINT_DBLTYPE const_6 = SIMINT_DBLSET1(6);
117     const SIMINT_DBLTYPE const_7 = SIMINT_DBLSET1(7);
118     const SIMINT_DBLTYPE const_8 = SIMINT_DBLSET1(8);
119     const SIMINT_DBLTYPE const_9 = SIMINT_DBLSET1(9);
120     const SIMINT_DBLTYPE one_half = SIMINT_DBLSET1(0.5);
121 
122 
123     ////////////////////////////////////////
124     // Loop over shells and primitives
125     ////////////////////////////////////////
126 
127     real_abcd = 0;
128     istart = 0;
129     for(ab = 0; ab < P.nshell12_clip; ++ab)
130     {
131         const int iend = istart + P.nprim12[ab];
132 
133         cd = 0;
134         jstart = 0;
135 
136         for(cd = 0; cd < Q.nshell12_clip; cd += SIMINT_NSHELL_SIMD)
137         {
138             const int nshellbatch = ((cd + SIMINT_NSHELL_SIMD) > Q.nshell12_clip) ? Q.nshell12_clip - cd : SIMINT_NSHELL_SIMD;
139             int jend = jstart;
140             for(i = 0; i < nshellbatch; i++)
141                 jend += Q.nprim12[cd+i];
142 
143             // Clear the beginning of the workspace (where we are accumulating integrals)
144             memset(work, 0, SIMINT_NSHELL_SIMD * 8280 * sizeof(double));
145             abcd = 0;
146 
147 
148             for(i = istart; i < iend; ++i)
149             {
150                 SIMINT_DBLTYPE bra_screen_max;  // only used if check_screen
151 
152                 if(check_screen)
153                 {
154                     // Skip this whole thing if always insignificant
155                     if((P.screen[i] * Q.screen_max) < screen_tol)
156                         continue;
157                     bra_screen_max = SIMINT_DBLSET1(P.screen[i]);
158                 }
159 
160                 icd = 0;
161                 iprimcd = 0;
162                 nprim_icd = Q.nprim12[cd];
163                 double * restrict PRIM_PTR_INT__i_s_g_s = INT__i_s_g_s + abcd * 420;
164                 double * restrict PRIM_PTR_INT__i_s_h_s = INT__i_s_h_s + abcd * 588;
165                 double * restrict PRIM_PTR_INT__k_s_g_s = INT__k_s_g_s + abcd * 540;
166                 double * restrict PRIM_PTR_INT__k_s_h_s = INT__k_s_h_s + abcd * 756;
167                 double * restrict PRIM_PTR_INT__l_s_g_s = INT__l_s_g_s + abcd * 675;
168                 double * restrict PRIM_PTR_INT__l_s_h_s = INT__l_s_h_s + abcd * 945;
169                 double * restrict PRIM_PTR_INT__m_s_g_s = INT__m_s_g_s + abcd * 825;
170                 double * restrict PRIM_PTR_INT__m_s_h_s = INT__m_s_h_s + abcd * 1155;
171                 double * restrict PRIM_PTR_INT__n_s_g_s = INT__n_s_g_s + abcd * 990;
172                 double * restrict PRIM_PTR_INT__n_s_h_s = INT__n_s_h_s + abcd * 1386;
173 
174 
175 
176                 // Load these one per loop over i
177                 const SIMINT_DBLTYPE P_alpha = SIMINT_DBLSET1(P.alpha[i]);
178                 const SIMINT_DBLTYPE P_prefac = SIMINT_DBLSET1(P.prefac[i]);
179                 const SIMINT_DBLTYPE Pxyz[3] = { SIMINT_DBLSET1(P.x[i]), SIMINT_DBLSET1(P.y[i]), SIMINT_DBLSET1(P.z[i]) };
180 
181                 const SIMINT_DBLTYPE P_PA[3] = { SIMINT_DBLSET1(P.PA_x[i]), SIMINT_DBLSET1(P.PA_y[i]), SIMINT_DBLSET1(P.PA_z[i]) };
182 
183                 for(j = jstart; j < jend; j += SIMINT_SIMD_LEN)
184                 {
185                     // calculate the shell offsets
186                     // these are the offset from the shell pointed to by cd
187                     // for each element
188                     int shelloffsets[SIMINT_SIMD_LEN] = {0};
189                     int lastoffset = 0;
190                     const int nlane = ( ((j + SIMINT_SIMD_LEN) < jend) ? SIMINT_SIMD_LEN : (jend - j));
191 
192                     if((iprimcd + SIMINT_SIMD_LEN) >= nprim_icd)
193                     {
194                         // Handle if the first element of the vector is a new shell
195                         if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
196                         {
197                             nprim_icd += Q.nprim12[cd + (++icd)];
198                             PRIM_PTR_INT__i_s_g_s += 420;
199                             PRIM_PTR_INT__i_s_h_s += 588;
200                             PRIM_PTR_INT__k_s_g_s += 540;
201                             PRIM_PTR_INT__k_s_h_s += 756;
202                             PRIM_PTR_INT__l_s_g_s += 675;
203                             PRIM_PTR_INT__l_s_h_s += 945;
204                             PRIM_PTR_INT__m_s_g_s += 825;
205                             PRIM_PTR_INT__m_s_h_s += 1155;
206                             PRIM_PTR_INT__n_s_g_s += 990;
207                             PRIM_PTR_INT__n_s_h_s += 1386;
208                         }
209                         iprimcd++;
210                         for(n = 1; n < SIMINT_SIMD_LEN; ++n)
211                         {
212                             if(iprimcd >= nprim_icd && ((icd+1) < nshellbatch))
213                             {
214                                 shelloffsets[n] = shelloffsets[n-1] + 1;
215                                 lastoffset++;
216                                 nprim_icd += Q.nprim12[cd + (++icd)];
217                             }
218                             else
219                                 shelloffsets[n] = shelloffsets[n-1];
220                             iprimcd++;
221                         }
222                     }
223                     else
224                         iprimcd += SIMINT_SIMD_LEN;
225 
226                     // Do we have to compute this vector (or has it been screened out)?
227                     // (not_screened != 0 means we have to do this vector)
228                     if(check_screen)
229                     {
230                         const double vmax = vector_max(SIMINT_MUL(bra_screen_max, SIMINT_DBLLOAD(Q.screen, j)));
231                         if(vmax < screen_tol)
232                         {
233                             PRIM_PTR_INT__i_s_g_s += lastoffset*420;
234                             PRIM_PTR_INT__i_s_h_s += lastoffset*588;
235                             PRIM_PTR_INT__k_s_g_s += lastoffset*540;
236                             PRIM_PTR_INT__k_s_h_s += lastoffset*756;
237                             PRIM_PTR_INT__l_s_g_s += lastoffset*675;
238                             PRIM_PTR_INT__l_s_h_s += lastoffset*945;
239                             PRIM_PTR_INT__m_s_g_s += lastoffset*825;
240                             PRIM_PTR_INT__m_s_h_s += lastoffset*1155;
241                             PRIM_PTR_INT__n_s_g_s += lastoffset*990;
242                             PRIM_PTR_INT__n_s_h_s += lastoffset*1386;
243                             continue;
244                         }
245                     }
246 
247                     const SIMINT_DBLTYPE Q_alpha = SIMINT_DBLLOAD(Q.alpha, j);
248                     const SIMINT_DBLTYPE PQalpha_mul = SIMINT_MUL(P_alpha, Q_alpha);
249                     const SIMINT_DBLTYPE PQalpha_sum = SIMINT_ADD(P_alpha, Q_alpha);
250                     const SIMINT_DBLTYPE one_over_PQalpha_sum = SIMINT_DIV(const_1, PQalpha_sum);
251 
252 
253                     /* construct R2 = (Px - Qx)**2 + (Py - Qy)**2 + (Pz -Qz)**2 */
254                     SIMINT_DBLTYPE PQ[3];
255                     PQ[0] = SIMINT_SUB(Pxyz[0], SIMINT_DBLLOAD(Q.x, j));
256                     PQ[1] = SIMINT_SUB(Pxyz[1], SIMINT_DBLLOAD(Q.y, j));
257                     PQ[2] = SIMINT_SUB(Pxyz[2], SIMINT_DBLLOAD(Q.z, j));
258                     SIMINT_DBLTYPE R2 = SIMINT_MUL(PQ[0], PQ[0]);
259                     R2 = SIMINT_FMADD(PQ[1], PQ[1], R2);
260                     R2 = SIMINT_FMADD(PQ[2], PQ[2], R2);
261 
262                     const SIMINT_DBLTYPE alpha = SIMINT_MUL(PQalpha_mul, one_over_PQalpha_sum); // alpha from MEST
263                     const SIMINT_DBLTYPE one_over_p = SIMINT_DIV(const_1, P_alpha);
264                     const SIMINT_DBLTYPE one_over_q = SIMINT_DIV(const_1, Q_alpha);
265                     const SIMINT_DBLTYPE one_over_2p = SIMINT_MUL(one_half, one_over_p);
266                     const SIMINT_DBLTYPE one_over_2q = SIMINT_MUL(one_half, one_over_q);
267                     const SIMINT_DBLTYPE one_over_2pq = SIMINT_MUL(one_half, one_over_PQalpha_sum);
268                     const SIMINT_DBLTYPE Q_PA[3] = { SIMINT_DBLLOAD(Q.PA_x, j), SIMINT_DBLLOAD(Q.PA_y, j), SIMINT_DBLLOAD(Q.PA_z, j) };
269 
270                     // NOTE: Minus sign!
271                     const SIMINT_DBLTYPE a_over_p = SIMINT_MUL(SIMINT_NEG(alpha), one_over_p);
272                     SIMINT_DBLTYPE aop_PQ[3];
273                     aop_PQ[0] = SIMINT_MUL(a_over_p, PQ[0]);
274                     aop_PQ[1] = SIMINT_MUL(a_over_p, PQ[1]);
275                     aop_PQ[2] = SIMINT_MUL(a_over_p, PQ[2]);
276 
277                     SIMINT_DBLTYPE a_over_q = SIMINT_MUL(alpha, one_over_q);
278                     SIMINT_DBLTYPE aoq_PQ[3];
279                     aoq_PQ[0] = SIMINT_MUL(a_over_q, PQ[0]);
280                     aoq_PQ[1] = SIMINT_MUL(a_over_q, PQ[1]);
281                     aoq_PQ[2] = SIMINT_MUL(a_over_q, PQ[2]);
282                     // Put a minus sign here so we don't have to in RR routines
283                     a_over_q = SIMINT_NEG(a_over_q);
284 
285 
286                     //////////////////////////////////////////////
287                     // Fjt function section
288                     // Maximum v value: 15
289                     //////////////////////////////////////////////
290                     // The parameter to the Fjt function
291                     const SIMINT_DBLTYPE F_x = SIMINT_MUL(R2, alpha);
292 
293 
294                     const SIMINT_DBLTYPE Q_prefac = mask_load(nlane, Q.prefac + j);
295 
296 
297                     boys_F_split(PRIM_INT__s_s_s_s, F_x, 15);
298                     SIMINT_DBLTYPE prefac = SIMINT_SQRT(one_over_PQalpha_sum);
299                     prefac = SIMINT_MUL(SIMINT_MUL(P_prefac, Q_prefac), prefac);
300                     for(n = 0; n <= 15; n++)
301                         PRIM_INT__s_s_s_s[n] = SIMINT_MUL(PRIM_INT__s_s_s_s[n], prefac);
302 
303                     //////////////////////////////////////////////
304                     // Primitive integrals: Vertical recurrance
305                     //////////////////////////////////////////////
306 
307                     const SIMINT_DBLTYPE vrr_const_1_over_2p = one_over_2p;
308                     const SIMINT_DBLTYPE vrr_const_2_over_2p = SIMINT_MUL(const_2, one_over_2p);
309                     const SIMINT_DBLTYPE vrr_const_3_over_2p = SIMINT_MUL(const_3, one_over_2p);
310                     const SIMINT_DBLTYPE vrr_const_4_over_2p = SIMINT_MUL(const_4, one_over_2p);
311                     const SIMINT_DBLTYPE vrr_const_5_over_2p = SIMINT_MUL(const_5, one_over_2p);
312                     const SIMINT_DBLTYPE vrr_const_6_over_2p = SIMINT_MUL(const_6, one_over_2p);
313                     const SIMINT_DBLTYPE vrr_const_7_over_2p = SIMINT_MUL(const_7, one_over_2p);
314                     const SIMINT_DBLTYPE vrr_const_8_over_2p = SIMINT_MUL(const_8, one_over_2p);
315                     const SIMINT_DBLTYPE vrr_const_9_over_2p = SIMINT_MUL(const_9, one_over_2p);
316                     const SIMINT_DBLTYPE vrr_const_1_over_2q = one_over_2q;
317                     const SIMINT_DBLTYPE vrr_const_2_over_2q = SIMINT_MUL(const_2, one_over_2q);
318                     const SIMINT_DBLTYPE vrr_const_3_over_2q = SIMINT_MUL(const_3, one_over_2q);
319                     const SIMINT_DBLTYPE vrr_const_4_over_2q = SIMINT_MUL(const_4, one_over_2q);
320                     const SIMINT_DBLTYPE vrr_const_1_over_2pq = one_over_2pq;
321                     const SIMINT_DBLTYPE vrr_const_2_over_2pq = SIMINT_MUL(const_2, one_over_2pq);
322                     const SIMINT_DBLTYPE vrr_const_3_over_2pq = SIMINT_MUL(const_3, one_over_2pq);
323                     const SIMINT_DBLTYPE vrr_const_4_over_2pq = SIMINT_MUL(const_4, one_over_2pq);
324                     const SIMINT_DBLTYPE vrr_const_5_over_2pq = SIMINT_MUL(const_5, one_over_2pq);
325                     const SIMINT_DBLTYPE vrr_const_6_over_2pq = SIMINT_MUL(const_6, one_over_2pq);
326                     const SIMINT_DBLTYPE vrr_const_7_over_2pq = SIMINT_MUL(const_7, one_over_2pq);
327                     const SIMINT_DBLTYPE vrr_const_8_over_2pq = SIMINT_MUL(const_8, one_over_2pq);
328                     const SIMINT_DBLTYPE vrr_const_9_over_2pq = SIMINT_MUL(const_9, one_over_2pq);
329                     const SIMINT_DBLTYPE vrr_const_10_over_2pq = SIMINT_MUL(const_10, one_over_2pq);
330 
331 
332 
333                     // Forming PRIM_INT__p_s_s_s[15 * 3];
334                     for(n = 0; n < 15; ++n)  // loop over orders of auxiliary function
335                     {
336 
337                         PRIM_INT__p_s_s_s[n * 3 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__s_s_s_s[n * 1 + 0]);
338                         PRIM_INT__p_s_s_s[n * 3 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__p_s_s_s[n * 3 + 0]);
339 
340                         PRIM_INT__p_s_s_s[n * 3 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__s_s_s_s[n * 1 + 0]);
341                         PRIM_INT__p_s_s_s[n * 3 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__p_s_s_s[n * 3 + 1]);
342 
343                         PRIM_INT__p_s_s_s[n * 3 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__s_s_s_s[n * 1 + 0]);
344                         PRIM_INT__p_s_s_s[n * 3 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__p_s_s_s[n * 3 + 2]);
345 
346                     }
347 
348 
349 
350                     // Forming PRIM_INT__d_s_s_s[14 * 6];
351                     for(n = 0; n < 14; ++n)  // loop over orders of auxiliary function
352                     {
353 
354                         PRIM_INT__d_s_s_s[n * 6 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__p_s_s_s[n * 3 + 0]);
355                         PRIM_INT__d_s_s_s[n * 6 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_s_s[n * 6 + 0]);
356                         PRIM_INT__d_s_s_s[n * 6 + 0] = SIMINT_FMADD( vrr_const_1_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__d_s_s_s[n * 6 + 0]);
357 
358                         PRIM_INT__d_s_s_s[n * 6 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__p_s_s_s[n * 3 + 0]);
359                         PRIM_INT__d_s_s_s[n * 6 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_s_s[n * 6 + 1]);
360 
361                         PRIM_INT__d_s_s_s[n * 6 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 0]);
362                         PRIM_INT__d_s_s_s[n * 6 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_s_s[n * 6 + 2]);
363 
364                         PRIM_INT__d_s_s_s[n * 6 + 3] = SIMINT_MUL(P_PA[1], PRIM_INT__p_s_s_s[n * 3 + 1]);
365                         PRIM_INT__d_s_s_s[n * 6 + 3] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_s_s[n * 6 + 3]);
366                         PRIM_INT__d_s_s_s[n * 6 + 3] = SIMINT_FMADD( vrr_const_1_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__d_s_s_s[n * 6 + 3]);
367 
368                         PRIM_INT__d_s_s_s[n * 6 + 4] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 1]);
369                         PRIM_INT__d_s_s_s[n * 6 + 4] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_s_s[n * 6 + 4]);
370 
371                         PRIM_INT__d_s_s_s[n * 6 + 5] = SIMINT_MUL(P_PA[2], PRIM_INT__p_s_s_s[n * 3 + 2]);
372                         PRIM_INT__d_s_s_s[n * 6 + 5] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__d_s_s_s[n * 6 + 5]);
373                         PRIM_INT__d_s_s_s[n * 6 + 5] = SIMINT_FMADD( vrr_const_1_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__s_s_s_s[(n+1) * 1 + 0], PRIM_INT__s_s_s_s[n * 1 + 0]), PRIM_INT__d_s_s_s[n * 6 + 5]);
374 
375                     }
376 
377 
378 
379                     // Forming PRIM_INT__f_s_s_s[13 * 10];
380                     for(n = 0; n < 13; ++n)  // loop over orders of auxiliary function
381                     {
382 
383                         PRIM_INT__f_s_s_s[n * 10 + 0] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 0]);
384                         PRIM_INT__f_s_s_s[n * 10 + 0] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__f_s_s_s[n * 10 + 0]);
385                         PRIM_INT__f_s_s_s[n * 10 + 0] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__p_s_s_s[n * 3 + 0]), PRIM_INT__f_s_s_s[n * 10 + 0]);
386 
387                         PRIM_INT__f_s_s_s[n * 10 + 1] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 0]);
388                         PRIM_INT__f_s_s_s[n * 10 + 1] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__f_s_s_s[n * 10 + 1]);
389 
390                         PRIM_INT__f_s_s_s[n * 10 + 2] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 0]);
391                         PRIM_INT__f_s_s_s[n * 10 + 2] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__f_s_s_s[n * 10 + 2]);
392 
393                         PRIM_INT__f_s_s_s[n * 10 + 3] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 3]);
394                         PRIM_INT__f_s_s_s[n * 10 + 3] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__f_s_s_s[n * 10 + 3]);
395 
396                         PRIM_INT__f_s_s_s[n * 10 + 4] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 1]);
397                         PRIM_INT__f_s_s_s[n * 10 + 4] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__f_s_s_s[n * 10 + 4]);
398 
399                         PRIM_INT__f_s_s_s[n * 10 + 5] = SIMINT_MUL(P_PA[0], PRIM_INT__d_s_s_s[n * 6 + 5]);
400                         PRIM_INT__f_s_s_s[n * 10 + 5] = SIMINT_FMADD( aop_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__f_s_s_s[n * 10 + 5]);
401 
402                         PRIM_INT__f_s_s_s[n * 10 + 6] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 3]);
403                         PRIM_INT__f_s_s_s[n * 10 + 6] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__f_s_s_s[n * 10 + 6]);
404                         PRIM_INT__f_s_s_s[n * 10 + 6] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__p_s_s_s[n * 3 + 1]), PRIM_INT__f_s_s_s[n * 10 + 6]);
405 
406                         PRIM_INT__f_s_s_s[n * 10 + 7] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 3]);
407                         PRIM_INT__f_s_s_s[n * 10 + 7] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__f_s_s_s[n * 10 + 7]);
408 
409                         PRIM_INT__f_s_s_s[n * 10 + 8] = SIMINT_MUL(P_PA[1], PRIM_INT__d_s_s_s[n * 6 + 5]);
410                         PRIM_INT__f_s_s_s[n * 10 + 8] = SIMINT_FMADD( aop_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__f_s_s_s[n * 10 + 8]);
411 
412                         PRIM_INT__f_s_s_s[n * 10 + 9] = SIMINT_MUL(P_PA[2], PRIM_INT__d_s_s_s[n * 6 + 5]);
413                         PRIM_INT__f_s_s_s[n * 10 + 9] = SIMINT_FMADD( aop_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__f_s_s_s[n * 10 + 9]);
414                         PRIM_INT__f_s_s_s[n * 10 + 9] = SIMINT_FMADD( vrr_const_2_over_2p, SIMINT_FMADD(a_over_p, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__p_s_s_s[n * 3 + 2]), PRIM_INT__f_s_s_s[n * 10 + 9]);
415 
416                     }
417 
418 
419                     VRR_I_g_s_s_s(
420                             PRIM_INT__g_s_s_s,
421                             PRIM_INT__f_s_s_s,
422                             PRIM_INT__d_s_s_s,
423                             P_PA,
424                             a_over_p,
425                             aop_PQ,
426                             one_over_2p,
427                             12);
428 
429 
430                     VRR_I_h_s_s_s(
431                             PRIM_INT__h_s_s_s,
432                             PRIM_INT__g_s_s_s,
433                             PRIM_INT__f_s_s_s,
434                             P_PA,
435                             a_over_p,
436                             aop_PQ,
437                             one_over_2p,
438                             11);
439 
440 
441                     ostei_general_vrr1_I(6, 10,
442                             one_over_2p, a_over_p, aop_PQ, P_PA,
443                             PRIM_INT__h_s_s_s, PRIM_INT__g_s_s_s, PRIM_INT__i_s_s_s);
444 
445 
446                     ostei_general_vrr_K(6, 0, 1, 0, 5,
447                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
448                             PRIM_INT__i_s_s_s, NULL, NULL, PRIM_INT__h_s_s_s, NULL, PRIM_INT__i_s_p_s);
449 
450 
451                     ostei_general_vrr_K(5, 0, 1, 0, 5,
452                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
453                             PRIM_INT__h_s_s_s, NULL, NULL, PRIM_INT__g_s_s_s, NULL, PRIM_INT__h_s_p_s);
454 
455 
456                     ostei_general_vrr_K(6, 0, 2, 0, 4,
457                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
458                             PRIM_INT__i_s_p_s, PRIM_INT__i_s_s_s, NULL, PRIM_INT__h_s_p_s, NULL, PRIM_INT__i_s_d_s);
459 
460 
461                     VRR_K_g_s_p_s(
462                             PRIM_INT__g_s_p_s,
463                             PRIM_INT__g_s_s_s,
464                             PRIM_INT__f_s_s_s,
465                             Q_PA,
466                             aoq_PQ,
467                             one_over_2pq,
468                             5);
469 
470 
471                     ostei_general_vrr_K(5, 0, 2, 0, 4,
472                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
473                             PRIM_INT__h_s_p_s, PRIM_INT__h_s_s_s, NULL, PRIM_INT__g_s_p_s, NULL, PRIM_INT__h_s_d_s);
474 
475 
476                     ostei_general_vrr_K(6, 0, 3, 0, 3,
477                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
478                             PRIM_INT__i_s_d_s, PRIM_INT__i_s_p_s, NULL, PRIM_INT__h_s_d_s, NULL, PRIM_INT__i_s_f_s);
479 
480 
481                     VRR_K_f_s_p_s(
482                             PRIM_INT__f_s_p_s,
483                             PRIM_INT__f_s_s_s,
484                             PRIM_INT__d_s_s_s,
485                             Q_PA,
486                             aoq_PQ,
487                             one_over_2pq,
488                             5);
489 
490 
491                     ostei_general_vrr_K(4, 0, 2, 0, 4,
492                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
493                             PRIM_INT__g_s_p_s, PRIM_INT__g_s_s_s, NULL, PRIM_INT__f_s_p_s, NULL, PRIM_INT__g_s_d_s);
494 
495 
496                     ostei_general_vrr_K(5, 0, 3, 0, 3,
497                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
498                             PRIM_INT__h_s_d_s, PRIM_INT__h_s_p_s, NULL, PRIM_INT__g_s_d_s, NULL, PRIM_INT__h_s_f_s);
499 
500 
501                     ostei_general_vrr_K(6, 0, 4, 0, 2,
502                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
503                             PRIM_INT__i_s_f_s, PRIM_INT__i_s_d_s, NULL, PRIM_INT__h_s_f_s, NULL, PRIM_INT__i_s_g_s);
504 
505 
506 
507                     // Forming PRIM_INT__d_s_p_s[5 * 18];
508                     for(n = 0; n < 5; ++n)  // loop over orders of auxiliary function
509                     {
510 
511                         PRIM_INT__d_s_p_s[n * 18 + 0] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 0]);
512                         PRIM_INT__d_s_p_s[n * 18 + 0] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__d_s_p_s[n * 18 + 0]);
513                         PRIM_INT__d_s_p_s[n * 18 + 0] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_p_s[n * 18 + 0]);
514 
515                         PRIM_INT__d_s_p_s[n * 18 + 1] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 0]);
516                         PRIM_INT__d_s_p_s[n * 18 + 1] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__d_s_p_s[n * 18 + 1]);
517 
518                         PRIM_INT__d_s_p_s[n * 18 + 2] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 0]);
519                         PRIM_INT__d_s_p_s[n * 18 + 2] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 0], PRIM_INT__d_s_p_s[n * 18 + 2]);
520 
521                         PRIM_INT__d_s_p_s[n * 18 + 3] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 1]);
522                         PRIM_INT__d_s_p_s[n * 18 + 3] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__d_s_p_s[n * 18 + 3]);
523                         PRIM_INT__d_s_p_s[n * 18 + 3] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_p_s[n * 18 + 3]);
524 
525                         PRIM_INT__d_s_p_s[n * 18 + 4] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 1]);
526                         PRIM_INT__d_s_p_s[n * 18 + 4] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__d_s_p_s[n * 18 + 4]);
527                         PRIM_INT__d_s_p_s[n * 18 + 4] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_p_s[n * 18 + 4]);
528 
529                         PRIM_INT__d_s_p_s[n * 18 + 5] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 1]);
530                         PRIM_INT__d_s_p_s[n * 18 + 5] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 1], PRIM_INT__d_s_p_s[n * 18 + 5]);
531 
532                         PRIM_INT__d_s_p_s[n * 18 + 6] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 2]);
533                         PRIM_INT__d_s_p_s[n * 18 + 6] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 2], PRIM_INT__d_s_p_s[n * 18 + 6]);
534                         PRIM_INT__d_s_p_s[n * 18 + 6] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__d_s_p_s[n * 18 + 6]);
535 
536                         PRIM_INT__d_s_p_s[n * 18 + 7] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 2]);
537                         PRIM_INT__d_s_p_s[n * 18 + 7] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 2], PRIM_INT__d_s_p_s[n * 18 + 7]);
538 
539                         PRIM_INT__d_s_p_s[n * 18 + 8] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 2]);
540                         PRIM_INT__d_s_p_s[n * 18 + 8] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 2], PRIM_INT__d_s_p_s[n * 18 + 8]);
541                         PRIM_INT__d_s_p_s[n * 18 + 8] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 0], PRIM_INT__d_s_p_s[n * 18 + 8]);
542 
543                         PRIM_INT__d_s_p_s[n * 18 + 9] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 3]);
544                         PRIM_INT__d_s_p_s[n * 18 + 9] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__d_s_p_s[n * 18 + 9]);
545 
546                         PRIM_INT__d_s_p_s[n * 18 + 10] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 3]);
547                         PRIM_INT__d_s_p_s[n * 18 + 10] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__d_s_p_s[n * 18 + 10]);
548                         PRIM_INT__d_s_p_s[n * 18 + 10] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_p_s[n * 18 + 10]);
549 
550                         PRIM_INT__d_s_p_s[n * 18 + 11] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 3]);
551                         PRIM_INT__d_s_p_s[n * 18 + 11] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 3], PRIM_INT__d_s_p_s[n * 18 + 11]);
552 
553                         PRIM_INT__d_s_p_s[n * 18 + 12] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 4]);
554                         PRIM_INT__d_s_p_s[n * 18 + 12] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 4], PRIM_INT__d_s_p_s[n * 18 + 12]);
555 
556                         PRIM_INT__d_s_p_s[n * 18 + 13] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 4]);
557                         PRIM_INT__d_s_p_s[n * 18 + 13] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 4], PRIM_INT__d_s_p_s[n * 18 + 13]);
558                         PRIM_INT__d_s_p_s[n * 18 + 13] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__d_s_p_s[n * 18 + 13]);
559 
560                         PRIM_INT__d_s_p_s[n * 18 + 14] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 4]);
561                         PRIM_INT__d_s_p_s[n * 18 + 14] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 4], PRIM_INT__d_s_p_s[n * 18 + 14]);
562                         PRIM_INT__d_s_p_s[n * 18 + 14] = SIMINT_FMADD( vrr_const_1_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 1], PRIM_INT__d_s_p_s[n * 18 + 14]);
563 
564                         PRIM_INT__d_s_p_s[n * 18 + 15] = SIMINT_MUL(Q_PA[0], PRIM_INT__d_s_s_s[n * 6 + 5]);
565                         PRIM_INT__d_s_p_s[n * 18 + 15] = SIMINT_FMADD( aoq_PQ[0], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__d_s_p_s[n * 18 + 15]);
566 
567                         PRIM_INT__d_s_p_s[n * 18 + 16] = SIMINT_MUL(Q_PA[1], PRIM_INT__d_s_s_s[n * 6 + 5]);
568                         PRIM_INT__d_s_p_s[n * 18 + 16] = SIMINT_FMADD( aoq_PQ[1], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__d_s_p_s[n * 18 + 16]);
569 
570                         PRIM_INT__d_s_p_s[n * 18 + 17] = SIMINT_MUL(Q_PA[2], PRIM_INT__d_s_s_s[n * 6 + 5]);
571                         PRIM_INT__d_s_p_s[n * 18 + 17] = SIMINT_FMADD( aoq_PQ[2], PRIM_INT__d_s_s_s[(n+1) * 6 + 5], PRIM_INT__d_s_p_s[n * 18 + 17]);
572                         PRIM_INT__d_s_p_s[n * 18 + 17] = SIMINT_FMADD( vrr_const_2_over_2pq, PRIM_INT__p_s_s_s[(n+1) * 3 + 2], PRIM_INT__d_s_p_s[n * 18 + 17]);
573 
574                     }
575 
576 
577                     VRR_K_f_s_d_s(
578                             PRIM_INT__f_s_d_s,
579                             PRIM_INT__f_s_p_s,
580                             PRIM_INT__f_s_s_s,
581                             PRIM_INT__d_s_p_s,
582                             Q_PA,
583                             a_over_q,
584                             aoq_PQ,
585                             one_over_2pq,
586                             one_over_2q,
587                             4);
588 
589 
590                     ostei_general_vrr_K(4, 0, 3, 0, 3,
591                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
592                             PRIM_INT__g_s_d_s, PRIM_INT__g_s_p_s, NULL, PRIM_INT__f_s_d_s, NULL, PRIM_INT__g_s_f_s);
593 
594 
595                     ostei_general_vrr_K(5, 0, 4, 0, 2,
596                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
597                             PRIM_INT__h_s_f_s, PRIM_INT__h_s_d_s, NULL, PRIM_INT__g_s_f_s, NULL, PRIM_INT__h_s_g_s);
598 
599 
600                     ostei_general_vrr_K(6, 0, 5, 0, 1,
601                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
602                             PRIM_INT__i_s_g_s, PRIM_INT__i_s_f_s, NULL, PRIM_INT__h_s_g_s, NULL, PRIM_INT__i_s_h_s);
603 
604 
605                     ostei_general_vrr1_I(7, 9,
606                             one_over_2p, a_over_p, aop_PQ, P_PA,
607                             PRIM_INT__i_s_s_s, PRIM_INT__h_s_s_s, PRIM_INT__k_s_s_s);
608 
609 
610                     ostei_general_vrr_K(7, 0, 1, 0, 5,
611                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
612                             PRIM_INT__k_s_s_s, NULL, NULL, PRIM_INT__i_s_s_s, NULL, PRIM_INT__k_s_p_s);
613 
614 
615                     ostei_general_vrr_K(7, 0, 2, 0, 4,
616                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
617                             PRIM_INT__k_s_p_s, PRIM_INT__k_s_s_s, NULL, PRIM_INT__i_s_p_s, NULL, PRIM_INT__k_s_d_s);
618 
619 
620                     ostei_general_vrr_K(7, 0, 3, 0, 3,
621                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
622                             PRIM_INT__k_s_d_s, PRIM_INT__k_s_p_s, NULL, PRIM_INT__i_s_d_s, NULL, PRIM_INT__k_s_f_s);
623 
624 
625                     ostei_general_vrr_K(7, 0, 4, 0, 2,
626                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
627                             PRIM_INT__k_s_f_s, PRIM_INT__k_s_d_s, NULL, PRIM_INT__i_s_f_s, NULL, PRIM_INT__k_s_g_s);
628 
629 
630                     ostei_general_vrr_K(7, 0, 5, 0, 1,
631                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
632                             PRIM_INT__k_s_g_s, PRIM_INT__k_s_f_s, NULL, PRIM_INT__i_s_g_s, NULL, PRIM_INT__k_s_h_s);
633 
634 
635                     ostei_general_vrr1_I(8, 8,
636                             one_over_2p, a_over_p, aop_PQ, P_PA,
637                             PRIM_INT__k_s_s_s, PRIM_INT__i_s_s_s, PRIM_INT__l_s_s_s);
638 
639 
640                     ostei_general_vrr_K(8, 0, 1, 0, 5,
641                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
642                             PRIM_INT__l_s_s_s, NULL, NULL, PRIM_INT__k_s_s_s, NULL, PRIM_INT__l_s_p_s);
643 
644 
645                     ostei_general_vrr_K(8, 0, 2, 0, 4,
646                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
647                             PRIM_INT__l_s_p_s, PRIM_INT__l_s_s_s, NULL, PRIM_INT__k_s_p_s, NULL, PRIM_INT__l_s_d_s);
648 
649 
650                     ostei_general_vrr_K(8, 0, 3, 0, 3,
651                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
652                             PRIM_INT__l_s_d_s, PRIM_INT__l_s_p_s, NULL, PRIM_INT__k_s_d_s, NULL, PRIM_INT__l_s_f_s);
653 
654 
655                     ostei_general_vrr_K(8, 0, 4, 0, 2,
656                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
657                             PRIM_INT__l_s_f_s, PRIM_INT__l_s_d_s, NULL, PRIM_INT__k_s_f_s, NULL, PRIM_INT__l_s_g_s);
658 
659 
660                     ostei_general_vrr_K(8, 0, 5, 0, 1,
661                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
662                             PRIM_INT__l_s_g_s, PRIM_INT__l_s_f_s, NULL, PRIM_INT__k_s_g_s, NULL, PRIM_INT__l_s_h_s);
663 
664 
665                     ostei_general_vrr1_I(9, 7,
666                             one_over_2p, a_over_p, aop_PQ, P_PA,
667                             PRIM_INT__l_s_s_s, PRIM_INT__k_s_s_s, PRIM_INT__m_s_s_s);
668 
669 
670                     ostei_general_vrr_K(9, 0, 1, 0, 5,
671                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
672                             PRIM_INT__m_s_s_s, NULL, NULL, PRIM_INT__l_s_s_s, NULL, PRIM_INT__m_s_p_s);
673 
674 
675                     ostei_general_vrr_K(9, 0, 2, 0, 4,
676                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
677                             PRIM_INT__m_s_p_s, PRIM_INT__m_s_s_s, NULL, PRIM_INT__l_s_p_s, NULL, PRIM_INT__m_s_d_s);
678 
679 
680                     ostei_general_vrr_K(9, 0, 3, 0, 3,
681                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
682                             PRIM_INT__m_s_d_s, PRIM_INT__m_s_p_s, NULL, PRIM_INT__l_s_d_s, NULL, PRIM_INT__m_s_f_s);
683 
684 
685                     ostei_general_vrr_K(9, 0, 4, 0, 2,
686                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
687                             PRIM_INT__m_s_f_s, PRIM_INT__m_s_d_s, NULL, PRIM_INT__l_s_f_s, NULL, PRIM_INT__m_s_g_s);
688 
689 
690                     ostei_general_vrr_K(9, 0, 5, 0, 1,
691                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
692                             PRIM_INT__m_s_g_s, PRIM_INT__m_s_f_s, NULL, PRIM_INT__l_s_g_s, NULL, PRIM_INT__m_s_h_s);
693 
694 
695                     ostei_general_vrr1_I(10, 6,
696                             one_over_2p, a_over_p, aop_PQ, P_PA,
697                             PRIM_INT__m_s_s_s, PRIM_INT__l_s_s_s, PRIM_INT__n_s_s_s);
698 
699 
700                     ostei_general_vrr_K(10, 0, 1, 0, 5,
701                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
702                             PRIM_INT__n_s_s_s, NULL, NULL, PRIM_INT__m_s_s_s, NULL, PRIM_INT__n_s_p_s);
703 
704 
705                     ostei_general_vrr_K(10, 0, 2, 0, 4,
706                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
707                             PRIM_INT__n_s_p_s, PRIM_INT__n_s_s_s, NULL, PRIM_INT__m_s_p_s, NULL, PRIM_INT__n_s_d_s);
708 
709 
710                     ostei_general_vrr_K(10, 0, 3, 0, 3,
711                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
712                             PRIM_INT__n_s_d_s, PRIM_INT__n_s_p_s, NULL, PRIM_INT__m_s_d_s, NULL, PRIM_INT__n_s_f_s);
713 
714 
715                     ostei_general_vrr_K(10, 0, 4, 0, 2,
716                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
717                             PRIM_INT__n_s_f_s, PRIM_INT__n_s_d_s, NULL, PRIM_INT__m_s_f_s, NULL, PRIM_INT__n_s_g_s);
718 
719 
720                     ostei_general_vrr_K(10, 0, 5, 0, 1,
721                             one_over_2q, a_over_q, one_over_2pq, aoq_PQ, Q_PA,
722                             PRIM_INT__n_s_g_s, PRIM_INT__n_s_f_s, NULL, PRIM_INT__m_s_g_s, NULL, PRIM_INT__n_s_h_s);
723 
724 
725 
726 
727                     ////////////////////////////////////
728                     // Accumulate contracted integrals
729                     ////////////////////////////////////
730                     if(lastoffset == 0)
731                     {
732                         contract_all(420, PRIM_INT__i_s_g_s, PRIM_PTR_INT__i_s_g_s);
733                         contract_all(588, PRIM_INT__i_s_h_s, PRIM_PTR_INT__i_s_h_s);
734                         contract_all(540, PRIM_INT__k_s_g_s, PRIM_PTR_INT__k_s_g_s);
735                         contract_all(756, PRIM_INT__k_s_h_s, PRIM_PTR_INT__k_s_h_s);
736                         contract_all(675, PRIM_INT__l_s_g_s, PRIM_PTR_INT__l_s_g_s);
737                         contract_all(945, PRIM_INT__l_s_h_s, PRIM_PTR_INT__l_s_h_s);
738                         contract_all(825, PRIM_INT__m_s_g_s, PRIM_PTR_INT__m_s_g_s);
739                         contract_all(1155, PRIM_INT__m_s_h_s, PRIM_PTR_INT__m_s_h_s);
740                         contract_all(990, PRIM_INT__n_s_g_s, PRIM_PTR_INT__n_s_g_s);
741                         contract_all(1386, PRIM_INT__n_s_h_s, PRIM_PTR_INT__n_s_h_s);
742                     }
743                     else
744                     {
745                         contract(420, shelloffsets, PRIM_INT__i_s_g_s, PRIM_PTR_INT__i_s_g_s);
746                         contract(588, shelloffsets, PRIM_INT__i_s_h_s, PRIM_PTR_INT__i_s_h_s);
747                         contract(540, shelloffsets, PRIM_INT__k_s_g_s, PRIM_PTR_INT__k_s_g_s);
748                         contract(756, shelloffsets, PRIM_INT__k_s_h_s, PRIM_PTR_INT__k_s_h_s);
749                         contract(675, shelloffsets, PRIM_INT__l_s_g_s, PRIM_PTR_INT__l_s_g_s);
750                         contract(945, shelloffsets, PRIM_INT__l_s_h_s, PRIM_PTR_INT__l_s_h_s);
751                         contract(825, shelloffsets, PRIM_INT__m_s_g_s, PRIM_PTR_INT__m_s_g_s);
752                         contract(1155, shelloffsets, PRIM_INT__m_s_h_s, PRIM_PTR_INT__m_s_h_s);
753                         contract(990, shelloffsets, PRIM_INT__n_s_g_s, PRIM_PTR_INT__n_s_g_s);
754                         contract(1386, shelloffsets, PRIM_INT__n_s_h_s, PRIM_PTR_INT__n_s_h_s);
755                         PRIM_PTR_INT__i_s_g_s += lastoffset*420;
756                         PRIM_PTR_INT__i_s_h_s += lastoffset*588;
757                         PRIM_PTR_INT__k_s_g_s += lastoffset*540;
758                         PRIM_PTR_INT__k_s_h_s += lastoffset*756;
759                         PRIM_PTR_INT__l_s_g_s += lastoffset*675;
760                         PRIM_PTR_INT__l_s_h_s += lastoffset*945;
761                         PRIM_PTR_INT__m_s_g_s += lastoffset*825;
762                         PRIM_PTR_INT__m_s_h_s += lastoffset*1155;
763                         PRIM_PTR_INT__n_s_g_s += lastoffset*990;
764                         PRIM_PTR_INT__n_s_h_s += lastoffset*1386;
765                     }
766 
767                 }  // close loop over j
768             }  // close loop over i
769 
770             //Advance to the next batch
771             jstart = SIMINT_SIMD_ROUND(jend);
772 
773             //////////////////////////////////////////////
774             // Contracted integrals: Horizontal recurrance
775             //////////////////////////////////////////////
776 
777 
778             const double hAB[3] = { P.AB_x[ab], P.AB_y[ab], P.AB_z[ab] };
779 
780 
781             for(abcd = 0; abcd < nshellbatch; ++abcd, ++real_abcd)
782             {
783                 const double hCD[3] = { Q.AB_x[cd+abcd], Q.AB_y[cd+abcd], Q.AB_z[cd+abcd] };
784 
785                 // set up HRR pointers
786                 double const * restrict HRR_INT__i_s_g_s = INT__i_s_g_s + abcd * 420;
787                 double const * restrict HRR_INT__i_s_h_s = INT__i_s_h_s + abcd * 588;
788                 double const * restrict HRR_INT__k_s_g_s = INT__k_s_g_s + abcd * 540;
789                 double const * restrict HRR_INT__k_s_h_s = INT__k_s_h_s + abcd * 756;
790                 double const * restrict HRR_INT__l_s_g_s = INT__l_s_g_s + abcd * 675;
791                 double const * restrict HRR_INT__l_s_h_s = INT__l_s_h_s + abcd * 945;
792                 double const * restrict HRR_INT__m_s_g_s = INT__m_s_g_s + abcd * 825;
793                 double const * restrict HRR_INT__m_s_h_s = INT__m_s_h_s + abcd * 1155;
794                 double const * restrict HRR_INT__n_s_g_s = INT__n_s_g_s + abcd * 990;
795                 double const * restrict HRR_INT__n_s_h_s = INT__n_s_h_s + abcd * 1386;
796                 double * restrict HRR_INT__i_g_g_p = INT__i_g_g_p + real_abcd * 18900;
797 
798                 // form INT__i_p_g_s
799                 ostei_general_hrr_J(6, 1, 4, 0, hAB, HRR_INT__k_s_g_s, HRR_INT__i_s_g_s, HRR_INT__i_p_g_s);
800 
801                 // form INT__i_p_h_s
802                 ostei_general_hrr_J(6, 1, 5, 0, hAB, HRR_INT__k_s_h_s, HRR_INT__i_s_h_s, HRR_INT__i_p_h_s);
803 
804                 // form INT__k_p_g_s
805                 ostei_general_hrr_J(7, 1, 4, 0, hAB, HRR_INT__l_s_g_s, HRR_INT__k_s_g_s, HRR_INT__k_p_g_s);
806 
807                 // form INT__k_p_h_s
808                 ostei_general_hrr_J(7, 1, 5, 0, hAB, HRR_INT__l_s_h_s, HRR_INT__k_s_h_s, HRR_INT__k_p_h_s);
809 
810                 // form INT__l_p_g_s
811                 ostei_general_hrr_J(8, 1, 4, 0, hAB, HRR_INT__m_s_g_s, HRR_INT__l_s_g_s, HRR_INT__l_p_g_s);
812 
813                 // form INT__l_p_h_s
814                 ostei_general_hrr_J(8, 1, 5, 0, hAB, HRR_INT__m_s_h_s, HRR_INT__l_s_h_s, HRR_INT__l_p_h_s);
815 
816                 // form INT__m_p_g_s
817                 ostei_general_hrr_J(9, 1, 4, 0, hAB, HRR_INT__n_s_g_s, HRR_INT__m_s_g_s, HRR_INT__m_p_g_s);
818 
819                 // form INT__m_p_h_s
820                 ostei_general_hrr_J(9, 1, 5, 0, hAB, HRR_INT__n_s_h_s, HRR_INT__m_s_h_s, HRR_INT__m_p_h_s);
821 
822                 // form INT__i_d_g_s
823                 ostei_general_hrr_J(6, 2, 4, 0, hAB, HRR_INT__k_p_g_s, HRR_INT__i_p_g_s, HRR_INT__i_d_g_s);
824 
825                 // form INT__i_d_h_s
826                 ostei_general_hrr_J(6, 2, 5, 0, hAB, HRR_INT__k_p_h_s, HRR_INT__i_p_h_s, HRR_INT__i_d_h_s);
827 
828                 // form INT__k_d_g_s
829                 ostei_general_hrr_J(7, 2, 4, 0, hAB, HRR_INT__l_p_g_s, HRR_INT__k_p_g_s, HRR_INT__k_d_g_s);
830 
831                 // form INT__k_d_h_s
832                 ostei_general_hrr_J(7, 2, 5, 0, hAB, HRR_INT__l_p_h_s, HRR_INT__k_p_h_s, HRR_INT__k_d_h_s);
833 
834                 // form INT__l_d_g_s
835                 ostei_general_hrr_J(8, 2, 4, 0, hAB, HRR_INT__m_p_g_s, HRR_INT__l_p_g_s, HRR_INT__l_d_g_s);
836 
837                 // form INT__l_d_h_s
838                 ostei_general_hrr_J(8, 2, 5, 0, hAB, HRR_INT__m_p_h_s, HRR_INT__l_p_h_s, HRR_INT__l_d_h_s);
839 
840                 // form INT__i_f_g_s
841                 ostei_general_hrr_J(6, 3, 4, 0, hAB, HRR_INT__k_d_g_s, HRR_INT__i_d_g_s, HRR_INT__i_f_g_s);
842 
843                 // form INT__i_f_h_s
844                 ostei_general_hrr_J(6, 3, 5, 0, hAB, HRR_INT__k_d_h_s, HRR_INT__i_d_h_s, HRR_INT__i_f_h_s);
845 
846                 // form INT__k_f_g_s
847                 ostei_general_hrr_J(7, 3, 4, 0, hAB, HRR_INT__l_d_g_s, HRR_INT__k_d_g_s, HRR_INT__k_f_g_s);
848 
849                 // form INT__k_f_h_s
850                 ostei_general_hrr_J(7, 3, 5, 0, hAB, HRR_INT__l_d_h_s, HRR_INT__k_d_h_s, HRR_INT__k_f_h_s);
851 
852                 // form INT__i_g_g_s
853                 ostei_general_hrr_J(6, 4, 4, 0, hAB, HRR_INT__k_f_g_s, HRR_INT__i_f_g_s, HRR_INT__i_g_g_s);
854 
855                 // form INT__i_g_h_s
856                 ostei_general_hrr_J(6, 4, 5, 0, hAB, HRR_INT__k_f_h_s, HRR_INT__i_f_h_s, HRR_INT__i_g_h_s);
857 
858                 // form INT__i_g_g_p
859                 HRR_L_g_p(
860                     HRR_INT__i_g_g_p,
861                     HRR_INT__i_g_g_s,
862                     HRR_INT__i_g_h_s,
863                     hCD, 420);
864 
865 
866             }  // close HRR loop
867 
868 
869         }   // close loop cdbatch
870 
871         istart = iend;
872     }  // close loop over ab
873 
874     return P.nshell12_clip * Q.nshell12_clip;
875 }
876 
ostei_g_i_g_p(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__g_i_g_p)877 int ostei_g_i_g_p(struct simint_multi_shellpair const P,
878                   struct simint_multi_shellpair const Q,
879                   double screen_tol,
880                   double * const restrict work,
881                   double * const restrict INT__g_i_g_p)
882 {
883     double P_AB[3*P.nshell12];
884     struct simint_multi_shellpair P_tmp = P;
885     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
886     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
887     P_tmp.AB_x = P_AB;
888     P_tmp.AB_y = P_AB + P.nshell12;
889     P_tmp.AB_z = P_AB + 2*P.nshell12;
890 
891     for(int i = 0; i < P.nshell12; i++)
892     {
893         P_tmp.AB_x[i] = -P.AB_x[i];
894         P_tmp.AB_y[i] = -P.AB_y[i];
895         P_tmp.AB_z[i] = -P.AB_z[i];
896     }
897 
898     int ret = ostei_i_g_g_p(P_tmp, Q, screen_tol, work, INT__g_i_g_p);
899     double buffer[18900] SIMINT_ALIGN_ARRAY_DBL;
900 
901     for(int q = 0; q < ret; q++)
902     {
903         int idx = 0;
904         for(int a = 0; a < 15; ++a)
905         for(int b = 0; b < 28; ++b)
906         for(int c = 0; c < 15; ++c)
907         for(int d = 0; d < 3; ++d)
908             buffer[idx++] = INT__g_i_g_p[q*18900+b*675+a*45+c*3+d];
909 
910         memcpy(INT__g_i_g_p+q*18900, buffer, 18900*sizeof(double));
911     }
912 
913     return ret;
914 }
915 
ostei_i_g_p_g(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__i_g_p_g)916 int ostei_i_g_p_g(struct simint_multi_shellpair const P,
917                   struct simint_multi_shellpair const Q,
918                   double screen_tol,
919                   double * const restrict work,
920                   double * const restrict INT__i_g_p_g)
921 {
922     double Q_AB[3*Q.nshell12];
923     struct simint_multi_shellpair Q_tmp = Q;
924     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
925     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
926     Q_tmp.AB_x = Q_AB;
927     Q_tmp.AB_y = Q_AB + Q.nshell12;
928     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
929 
930     for(int i = 0; i < Q.nshell12; i++)
931     {
932         Q_tmp.AB_x[i] = -Q.AB_x[i];
933         Q_tmp.AB_y[i] = -Q.AB_y[i];
934         Q_tmp.AB_z[i] = -Q.AB_z[i];
935     }
936 
937     int ret = ostei_i_g_g_p(P, Q_tmp, screen_tol, work, INT__i_g_p_g);
938     double buffer[18900] SIMINT_ALIGN_ARRAY_DBL;
939 
940     for(int q = 0; q < ret; q++)
941     {
942         int idx = 0;
943         for(int a = 0; a < 28; ++a)
944         for(int b = 0; b < 15; ++b)
945         for(int c = 0; c < 3; ++c)
946         for(int d = 0; d < 15; ++d)
947             buffer[idx++] = INT__i_g_p_g[q*18900+a*675+b*45+d*3+c];
948 
949         memcpy(INT__i_g_p_g+q*18900, buffer, 18900*sizeof(double));
950     }
951 
952     return ret;
953 }
954 
ostei_g_i_p_g(struct simint_multi_shellpair const P,struct simint_multi_shellpair const Q,double screen_tol,double * const restrict work,double * const restrict INT__g_i_p_g)955 int ostei_g_i_p_g(struct simint_multi_shellpair const P,
956                   struct simint_multi_shellpair const Q,
957                   double screen_tol,
958                   double * const restrict work,
959                   double * const restrict INT__g_i_p_g)
960 {
961     double P_AB[3*P.nshell12];
962     struct simint_multi_shellpair P_tmp = P;
963     P_tmp.PA_x = P.PB_x;  P_tmp.PA_y = P.PB_y;  P_tmp.PA_z = P.PB_z;
964     P_tmp.PB_x = P.PA_x;  P_tmp.PB_y = P.PA_y;  P_tmp.PB_z = P.PA_z;
965     P_tmp.AB_x = P_AB;
966     P_tmp.AB_y = P_AB + P.nshell12;
967     P_tmp.AB_z = P_AB + 2*P.nshell12;
968 
969     for(int i = 0; i < P.nshell12; i++)
970     {
971         P_tmp.AB_x[i] = -P.AB_x[i];
972         P_tmp.AB_y[i] = -P.AB_y[i];
973         P_tmp.AB_z[i] = -P.AB_z[i];
974     }
975 
976     double Q_AB[3*Q.nshell12];
977     struct simint_multi_shellpair Q_tmp = Q;
978     Q_tmp.PA_x = Q.PB_x;  Q_tmp.PA_y = Q.PB_y;  Q_tmp.PA_z = Q.PB_z;
979     Q_tmp.PB_x = Q.PA_x;  Q_tmp.PB_y = Q.PA_y;  Q_tmp.PB_z = Q.PA_z;
980     Q_tmp.AB_x = Q_AB;
981     Q_tmp.AB_y = Q_AB + Q.nshell12;
982     Q_tmp.AB_z = Q_AB + 2*Q.nshell12;
983 
984     for(int i = 0; i < Q.nshell12; i++)
985     {
986         Q_tmp.AB_x[i] = -Q.AB_x[i];
987         Q_tmp.AB_y[i] = -Q.AB_y[i];
988         Q_tmp.AB_z[i] = -Q.AB_z[i];
989     }
990 
991     int ret = ostei_i_g_g_p(P_tmp, Q_tmp, screen_tol, work, INT__g_i_p_g);
992     double buffer[18900] SIMINT_ALIGN_ARRAY_DBL;
993 
994     for(int q = 0; q < ret; q++)
995     {
996         int idx = 0;
997         for(int a = 0; a < 15; ++a)
998         for(int b = 0; b < 28; ++b)
999         for(int c = 0; c < 3; ++c)
1000         for(int d = 0; d < 15; ++d)
1001             buffer[idx++] = INT__g_i_p_g[q*18900+b*675+a*45+d*3+c];
1002 
1003         memcpy(INT__g_i_p_g+q*18900, buffer, 18900*sizeof(double));
1004     }
1005 
1006     return ret;
1007 }
1008 
1009