src/mdlib/nbnxn_kernels/nbnxn_kernel_simd_utils.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
   5  * Copyright (c) 2001-2012, The GROMACS Development Team
   6  * Copyright (c) 2012,2013, by the GROMACS development team, led by
   7  * David van der Spoel, Berk Hess, Erik Lindahl, and including many
   8  * others, as listed in the AUTHORS file in the top-level source
   9  * directory and at http://www.gromacs.org.
  10  *
  11  * GROMACS is free software; you can redistribute it and/or
  12  * modify it under the terms of the GNU Lesser General Public License
  13  * as published by the Free Software Foundation; either version 2.1
  14  * of the License, or (at your option) any later version.
  15  *
  16  * GROMACS is distributed in the hope that it will be useful,
  17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19  * Lesser General Public License for more details.
  20  *
  21  * You should have received a copy of the GNU Lesser General Public
  22  * License along with GROMACS; if not, see
  23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  25  *
  26  * If you want to redistribute modifications to GROMACS, please
  27  * consider that scientific software is very special. Version
  28  * control is crucial - bugs must be traceable. We will be happy to
  29  * consider code for inclusion in the official distribution, but
  30  * derived work must not be called official GROMACS. Details are found
  31  * in the README & COPYING files - if they are missing, get the
  32  * official version at http://www.gromacs.org.
  33  *
  34  * To help us fund GROMACS development, we humbly ask that you cite
  35  * the research papers on the package. Check out http://www.gromacs.org.
  36  */
  37 #ifndef _nbnxn_kernel_sse_utils_h_
  38 #define _nbnxn_kernel_sse_utils_h_
  39
  40 /* This files contains all functions/macros for the SIMD kernels
  41  * which have explicit dependencies on the j-cluster size and/or SIMD-width.
  42  * The functionality which depends on the j-cluster size is:
  43  *   LJ-parameter lookup
  44  *   force table lookup
  45  *   energy group pair energy storage
  46  */
  47
  48 #ifdef GMX_X86_SSE2
  49
  50 /* Transpose 2 double precision registers */
  51 #define GMX_MM_TRANSPOSE2_OP_PD(in0,in1,out0,out1)                      \
  52 {                                                                       \
  53     out0 = _mm_unpacklo_pd(in0,in1);                                    \
  54     out1 = _mm_unpackhi_pd(in0,in1);                                    \
  55 }
  56
  57 #if defined GMX_MM128_HERE || !defined GMX_DOUBLE
  58 /* Collect element 0 and 1 of the 4 inputs to out0 and out1, respectively */
  59 #define GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(in0,in1,in2,in3,out0,out1)    \
  60 {                                                                       \
  61     __m128 _c01,_c23;                                                   \
  62     _c01 = _mm_movelh_ps(in0,in1);                                      \
  63     _c23 = _mm_movelh_ps(in2,in3);                                      \
  64     out0 = _mm_shuffle_ps(_c01,_c23,_MM_SHUFFLE(2,0,2,0));              \
  65     out1 = _mm_shuffle_ps(_c01,_c23,_MM_SHUFFLE(3,1,3,1));              \
  66 }
  67 #else
  68 /* Collect element 0 and 1 of the 4 inputs to out0 and out1, respectively */
  69 #define GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(in0,in1,in2,in3,out0,out1)    \
  70 {                                                                       \
  71     __m256d _c01,_c23;                                                  \
  72     _c01 = _mm256_shuffle_pd(in0,in1,_MM_SHUFFLE(1,0,1,0));             \
  73     _c23 = _mm256_shuffle_pd(in2,in3,_MM_SHUFFLE(1,0,1,0));             \
  74     out0 = _mm256_shuffle_pd(_c01,_c23,_MM_SHUFFLE(2,0,2,0));           \
  75     out1 = _mm256_shuffle_pd(_c01,_c23,_MM_SHUFFLE(3,1,3,1));           \
  76 }
  77 #endif
  78
  79 /* Collect element 2 of the 4 inputs to out */
  80 #define GMX_MM_SHUFFLE_4_PS_FIL2_TO_1_PS(in0,in1,in2,in3,out)           \
  81 {                                                                       \
  82     __m128 _c01,_c23;                                                   \
  83     _c01 = _mm_shuffle_ps(in0,in1,_MM_SHUFFLE(3,2,3,2));                \
  84     _c23 = _mm_shuffle_ps(in2,in3,_MM_SHUFFLE(3,2,3,2));                \
  85     out  = _mm_shuffle_ps(_c01,_c23,_MM_SHUFFLE(2,0,2,0));              \
  86 }
  87
  88 #ifndef GMX_MM256_HERE
  89 #ifndef GMX_DOUBLE
  90 /* Sum the elements within each input register and store the sums in out */
  91 #define GMX_MM_TRANSPOSE_SUM4_PR(in0,in1,in2,in3,out)                   \
  92 {                                                                       \
  93     _MM_TRANSPOSE4_PS(in0,in1,in2,in3);                                 \
  94     in0 = _mm_add_ps(in0,in1);                                          \
  95     in2 = _mm_add_ps(in2,in3);                                          \
  96     out  = _mm_add_ps(in0,in2);                                         \
  97 }
  98 #else
  99 /* Sum the elements within each input register and store the sums in out */
 100 #define GMX_MM_TRANSPOSE_SUM2_PD(in0,in1,out)                           \
 101 {                                                                       \
 102     GMX_MM_TRANSPOSE2_PD(in0,in1);                                      \
 103     out  = _mm_add_pd(in0,in1);                                         \
 104 }
 105 #endif
 106 #else
 107 #ifndef GMX_DOUBLE
 108 /* Sum the elements within each input register and store the sums in out */
 109 #define GMX_MM_TRANSPOSE_SUM4_PR(in0,in1,in2,in3,out)                   \
 110 {                                                                       \
 111     in0 = _mm256_hadd_ps(in0,in1);                                      \
 112     in2 = _mm256_hadd_ps(in2,in3);                                      \
 113     in1 = _mm256_hadd_ps(in0,in2);                                      \
 114     out = _mm_add_ps(_mm256_castps256_ps128(in1),_mm256_extractf128_ps(in1,1)); \
 115 }
 116 /* Sum the elements of halfs of each input register and store sums in out */
 117 #define GMX_MM_TRANSPOSE_SUM4H_PR(in0,in2,out)                          \
 118 {                                                                       \
 119     in0 = _mm256_hadd_ps(in0,_mm256_setzero_ps());                      \
 120     in2 = _mm256_hadd_ps(in2,_mm256_setzero_ps());                      \
 121     in0 = _mm256_hadd_ps(in0,in2);                                      \
 122     in2 = _mm256_permute_ps(in0,_MM_SHUFFLE(2,3,0,1));                  \
 123     out = _mm_add_ps(_mm256_castps256_ps128(in0),_mm256_extractf128_ps(in2,1)); \
 124 }
 125 #else
 126 /* Sum the elements within each input register and store the sums in out */
 127 #define GMX_MM_TRANSPOSE_SUM4_PR(in0,in1,in2,in3,out)                   \
 128 {                                                                       \
 129     in0 = _mm256_hadd_pd(in0,in1);                                      \
 130     in2 = _mm256_hadd_pd(in2,in3);                                      \
 131     out = _mm256_add_pd(_mm256_permute2f128_pd(in0,in2,0x20),_mm256_permute2f128_pd(in0,in2,0x31)); \
 132 }
 133 #endif
 134 #endif
 135
 136 #ifdef GMX_MM128_HERE
 137
 138 static inline __m128
 139 gmx_mm128_invsqrt_ps_single(__m128 x)
 140 {
 141     const __m128 half  = _mm_set_ps(0.5,0.5,0.5,0.5);
 142     const __m128 three = _mm_set_ps(3.0,3.0,3.0,3.0);
 143
 144     __m128 lu = _mm_rsqrt_ps(x);
 145
 146     return _mm_mul_ps(half,_mm_mul_ps(_mm_sub_ps(three,_mm_mul_ps(_mm_mul_ps(lu,lu),x)),lu));
 147 }
 148
 149 /* Do 2 double precision invsqrt operations.
 150  * Doing the SIMD rsqrt and the first Newton Raphson iteration
 151  * in single precision gives full double precision accuracy.
 152  * The speed is more than double that of two gmx_mm_invsqrt_pd calls.
 153  */
 154 #define GMX_MM128_INVSQRT2_PD(in0,in1,out0,out1)                        \
 155 {                                                                       \
 156     const __m128d half  = _mm_set1_pd(0.5);                             \
 157     const __m128d three = _mm_set1_pd(3.0);                             \
 158     __m128  s,ir;                                                       \
 159     __m128d lu0,lu1;                                                    \
 160                                                                         \
 161     s    = _mm_movelh_ps(_mm_cvtpd_ps(in0),_mm_cvtpd_ps(in1));          \
 162     ir   = gmx_mm128_invsqrt_ps_single(s);                              \
 163     lu0  = _mm_cvtps_pd(ir);                                            \
 164     lu1  = _mm_cvtps_pd(_mm_movehl_ps(ir,ir));                          \
 165     out0 = _mm_mul_pd(half,_mm_mul_pd(_mm_sub_pd(three,_mm_mul_pd(_mm_mul_pd(lu0,lu0),in0)),lu0)); \
 166     out1 = _mm_mul_pd(half,_mm_mul_pd(_mm_sub_pd(three,_mm_mul_pd(_mm_mul_pd(lu1,lu1),in1)),lu1)); \
 167 }
 168
 169 #define GMX_MM_INVSQRT2_PD GMX_MM128_INVSQRT2_PD
 170
 171 #endif
 172
 173 #ifdef GMX_MM256_HERE
 174
 175 static inline __m256
 176 gmx_mm256_invsqrt_ps_single(__m256 x)
 177 {
 178     const __m256 half  = _mm256_set_ps(0.5,0.5,0.5,0.5,0.5,0.5,0.5,0.5);
 179     const __m256 three = _mm256_set_ps(3.0,3.0,3.0,3.0,3.0,3.0,3.0,3.0);
 180
 181     __m256 lu = _mm256_rsqrt_ps(x);
 182
 183     return _mm256_mul_ps(half,_mm256_mul_ps(_mm256_sub_ps(three,_mm256_mul_ps(_mm256_mul_ps(lu,lu),x)),lu));
 184 }
 185
 186 /* Do 4 double precision invsqrt operations.
 187  * Doing the SIMD rsqrt and the first Newton Raphson iteration
 188  * in single precision gives full double precision accuracy.
 189  */
 190 #define GMX_MM256_INVSQRT2_PD(in0,in1,out0,out1)                        \
 191 {                                                                       \
 192     const __m256d half  = _mm256_set1_pd(0.5);                          \
 193     const __m256d three = _mm256_set1_pd(3.0);                          \
 194     __m256  s,ir;                                                       \
 195     __m256d lu0,lu1;                                                    \
 196                                                                         \
 197     s    = _mm256_insertf128_ps(_mm256_castps128_ps256(_mm256_cvtpd_ps(in0)),_mm256_cvtpd_ps(in1),1); \
 198     ir   = gmx_mm256_invsqrt_ps_single(s);                              \
 199     lu0  = _mm256_cvtps_pd(_mm256_castps256_ps128(ir));                 \
 200     lu1  = _mm256_cvtps_pd(_mm256_extractf128_ps(ir,1));                \
 201     out0 = _mm256_mul_pd(half,_mm256_mul_pd(_mm256_sub_pd(three,_mm256_mul_pd(_mm256_mul_pd(lu0,lu0),in0)),lu0)); \
 202     out1 = _mm256_mul_pd(half,_mm256_mul_pd(_mm256_sub_pd(three,_mm256_mul_pd(_mm256_mul_pd(lu1,lu1),in1)),lu1)); \
 203 }
 204
 205 #define GMX_MM_INVSQRT2_PD GMX_MM256_INVSQRT2_PD
 206
 207 #endif
 208
 209 /* Force and energy table load and interpolation routines */
 210
 211 #if defined GMX_MM128_HERE && !defined GMX_DOUBLE
 212
 213 #define load_lj_pair_params(nbfp,type,aj,c6_SSE,c12_SSE)                \
 214 {                                                                       \
 215     gmx_mm_pr clj_SSE[UNROLLJ];                                         \
 216     int p;                                                              \
 217                                                                         \
 218     for(p=0; p<UNROLLJ; p++)                                            \
 219     {                                                                   \
 220         /* Here we load 4 aligned floats, but we need just 2 */         \
 221         clj_SSE[p] = gmx_load_pr(nbfp+type[aj+p]*NBFP_STRIDE);          \
 222     }                                                                   \
 223     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(clj_SSE[0],clj_SSE[1],clj_SSE[2],clj_SSE[3],c6_SSE,c12_SSE); \
 224 }
 225
 226 #endif
 227
 228 #if defined GMX_MM256_HERE && !defined GMX_DOUBLE
 229
 230 /* Put two 128-bit 4-float registers into one 256-bit 8-float register */
 231 #define GMX_2_MM_TO_M256(in0,in1,out)                                   \
 232 {                                                                       \
 233     out = _mm256_insertf128_ps(_mm256_castps128_ps256(in0),in1,1);      \
 234 }
 235
 236 #define load_lj_pair_params(nbfp,type,aj,c6_SSE,c12_SSE)                \
 237 {                                                                       \
 238     __m128 clj_SSE[UNROLLJ],c6t_SSE[2],c12t_SSE[2];                     \
 239     int p;                                                              \
 240                                                                         \
 241     for(p=0; p<UNROLLJ; p++)                                            \
 242     {                                                                   \
 243         /* Here we load 4 aligned floats, but we need just 2 */         \
 244         clj_SSE[p] = _mm_load_ps(nbfp+type[aj+p]*NBFP_STRIDE);          \
 245     }                                                                   \
 246     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(clj_SSE[0],clj_SSE[1],clj_SSE[2],clj_SSE[3],c6t_SSE[0],c12t_SSE[0]); \
 247     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(clj_SSE[4],clj_SSE[5],clj_SSE[6],clj_SSE[7],c6t_SSE[1],c12t_SSE[1]); \
 248                                                                         \
 249     GMX_2_MM_TO_M256(c6t_SSE[0],c6t_SSE[1],c6_SSE);                     \
 250     GMX_2_MM_TO_M256(c12t_SSE[0],c12t_SSE[1],c12_SSE);                  \
 251 }
 252
 253 #define load_lj_pair_params2(nbfp0,nbfp1,type,aj,c6_SSE,c12_SSE)        \
 254 {                                                                       \
 255     __m128 clj_SSE0[UNROLLJ],clj_SSE1[UNROLLJ],c6t_SSE[2],c12t_SSE[2];  \
 256     int p;                                                              \
 257                                                                         \
 258     for(p=0; p<UNROLLJ; p++)                                            \
 259     {                                                                   \
 260         /* Here we load 4 aligned floats, but we need just 2 */         \
 261         clj_SSE0[p] = _mm_load_ps(nbfp0+type[aj+p]*NBFP_STRIDE);        \
 262     }                                                                   \
 263     for(p=0; p<UNROLLJ; p++)                                            \
 264     {                                                                   \
 265         /* Here we load 4 aligned floats, but we need just 2 */         \
 266         clj_SSE1[p] = _mm_load_ps(nbfp1+type[aj+p]*NBFP_STRIDE);        \
 267     }                                                                   \
 268     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(clj_SSE0[0],clj_SSE0[1],clj_SSE0[2],clj_SSE0[3],c6t_SSE[0],c12t_SSE[0]); \
 269     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(clj_SSE1[0],clj_SSE1[1],clj_SSE1[2],clj_SSE1[3],c6t_SSE[1],c12t_SSE[1]); \
 270                                                                         \
 271     GMX_2_MM_TO_M256(c6t_SSE[0],c6t_SSE[1],c6_SSE);                     \
 272     GMX_2_MM_TO_M256(c12t_SSE[0],c12t_SSE[1],c12_SSE);                  \
 273 }
 274
 275 #endif
 276
 277 #if defined GMX_MM128_HERE && defined GMX_DOUBLE
 278
 279 #define load_lj_pair_params(nbfp,type,aj,c6_SSE,c12_SSE)                \
 280 {                                                                       \
 281     gmx_mm_pr clj_SSE[UNROLLJ];                                         \
 282     int p;                                                              \
 283                                                                         \
 284     for(p=0; p<UNROLLJ; p++)                                            \
 285     {                                                                   \
 286         clj_SSE[p] = gmx_load_pr(nbfp+type[aj+p]*NBFP_STRIDE);          \
 287     }                                                                   \
 288     GMX_MM_TRANSPOSE2_OP_PD(clj_SSE[0],clj_SSE[1],c6_SSE,c12_SSE);      \
 289 }
 290
 291 #endif
 292
 293 #if defined GMX_MM256_HERE && defined GMX_DOUBLE
 294
 295 #define load_lj_pair_params(nbfp,type,aj,c6_SSE,c12_SSE)                \
 296 {                                                                       \
 297     __m128d clj_SSE[UNROLLJ],c6t_SSE[2],c12t_SSE[2];                    \
 298     int p;                                                              \
 299                                                                         \
 300     for(p=0; p<UNROLLJ; p++)                                            \
 301     {                                                                   \
 302         clj_SSE[p] = _mm_load_pd(nbfp+type[aj+p]*NBFP_STRIDE);          \
 303     }                                                                   \
 304     GMX_MM_TRANSPOSE2_OP_PD(clj_SSE[0],clj_SSE[1],c6t_SSE[0],c12t_SSE[0]); \
 305     GMX_MM_TRANSPOSE2_OP_PD(clj_SSE[2],clj_SSE[3],c6t_SSE[1],c12t_SSE[1]); \
 306     GMX_2_M128D_TO_M256D(c6t_SSE[0],c6t_SSE[1],c6_SSE);                 \
 307     GMX_2_M128D_TO_M256D(c12t_SSE[0],c12t_SSE[1],c12_SSE);              \
 308 }
 309
 310 #endif
 311
 312
 313 /* The load_table functions below are performance critical.
 314  * The routines issue UNROLLI*UNROLLJ _mm_load_ps calls.
 315  * As these all have latencies, scheduling is crucial.
 316  * The Intel compilers and CPUs seem to do a good job at this.
 317  * But AMD CPUs perform significantly worse with gcc than with icc.
 318  * Performance is improved a bit by using the extract function UNROLLJ times,
 319  * instead of doing an _mm_store_si128 for every i-particle.
 320  * This is only faster when we use FDV0 formatted tables, where we also need
 321  * to multiple the index by 4, which can be done by a SIMD bit shift.
 322  * With single precision AVX, 8 extracts are much slower than 1 store.
 323  * Because of this, the load_table_f macro always takes the ti parameter,
 324  * but it is only used with AVX.
 325  */
 326
 327 #if defined GMX_MM128_HERE && !defined GMX_DOUBLE
 328
 329 #define load_table_f(tab_coul_FDV0, ti_SSE, ti, ctab0_SSE, ctab1_SSE)   \
 330 {                                                                       \
 331     int    idx[4];                                                      \
 332     __m128 ctab_SSE[4];                                                 \
 333                                                                         \
 334     /* Table has 4 entries, left-shift index by 2 */                    \
 335     ti_SSE = _mm_slli_epi32(ti_SSE,2);                                  \
 336     /* Without SSE4.1 the extract macro needs an immediate: unroll */   \
 337     idx[0] = gmx_mm_extract_epi32(ti_SSE,0);                            \
 338     ctab_SSE[0] = _mm_load_ps(tab_coul_FDV0+idx[0]);                    \
 339     idx[1] = gmx_mm_extract_epi32(ti_SSE,1);                            \
 340     ctab_SSE[1] = _mm_load_ps(tab_coul_FDV0+idx[1]);                    \
 341     idx[2] = gmx_mm_extract_epi32(ti_SSE,2);                            \
 342     ctab_SSE[2] = _mm_load_ps(tab_coul_FDV0+idx[2]);                    \
 343     idx[3] = gmx_mm_extract_epi32(ti_SSE,3);                            \
 344     ctab_SSE[3] = _mm_load_ps(tab_coul_FDV0+idx[3]);                    \
 345                                                                         \
 346     /* Shuffle the force table entries to a convenient order */         \
 347     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(ctab_SSE[0],ctab_SSE[1],ctab_SSE[2],ctab_SSE[3],ctab0_SSE,ctab1_SSE); \
 348 }
 349
 350 #define load_table_f_v(tab_coul_FDV0, ti_SSE, ti, ctab0_SSE, ctab1_SSE, ctabv_SSE) \
 351 {                                                                       \
 352     int    idx[4];                                                      \
 353     __m128 ctab_SSE[4];                                                 \
 354                                                                         \
 355     /* Table has 4 entries, left-shift index by 2 */                    \
 356     ti_SSE = _mm_slli_epi32(ti_SSE,2);                                  \
 357     /* Without SSE4.1 the extract macro needs an immediate: unroll */   \
 358     idx[0] = gmx_mm_extract_epi32(ti_SSE,0);                            \
 359     ctab_SSE[0] = _mm_load_ps(tab_coul_FDV0+idx[0]);                    \
 360     idx[1] = gmx_mm_extract_epi32(ti_SSE,1);                            \
 361     ctab_SSE[1] = _mm_load_ps(tab_coul_FDV0+idx[1]);                    \
 362     idx[2] = gmx_mm_extract_epi32(ti_SSE,2);                            \
 363     ctab_SSE[2] = _mm_load_ps(tab_coul_FDV0+idx[2]);                    \
 364     idx[3] = gmx_mm_extract_epi32(ti_SSE,3);                            \
 365     ctab_SSE[3] = _mm_load_ps(tab_coul_FDV0+idx[3]);                    \
 366                                                                         \
 367     /* Shuffle the force  table entries to a convenient order */        \
 368     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(ctab_SSE[0],ctab_SSE[1],ctab_SSE[2],ctab_SSE[3],ctab0_SSE,ctab1_SSE); \
 369     /* Shuffle the energy table entries to a convenient order */        \
 370     GMX_MM_SHUFFLE_4_PS_FIL2_TO_1_PS(ctab_SSE[0],ctab_SSE[1],ctab_SSE[2],ctab_SSE[3],ctabv_SSE); \
 371 }
 372
 373 #endif
 374
 375 #if defined GMX_MM256_HERE && !defined GMX_DOUBLE
 376
 377 #define load_table_f(tab_coul_FDV0, ti_SSE, ti, ctab0_SSE, ctab1_SSE)   \
 378 {                                                                       \
 379     __m128 ctab_SSE[8],ctabt_SSE[4];                                    \
 380     int    j;                                                           \
 381                                                                         \
 382     /* Bit shifting would be faster, but AVX doesn't support that */    \
 383     _mm256_store_si256((__m256i *)ti,ti_SSE);                           \
 384     for(j=0; j<8; j++)                                                  \
 385     {                                                                   \
 386         ctab_SSE[j] = _mm_load_ps(tab_coul_FDV0+ti[j]*4);               \
 387     }                                                                   \
 388     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(ctab_SSE[0],ctab_SSE[1],ctab_SSE[2],ctab_SSE[3],ctabt_SSE[0],ctabt_SSE[2]); \
 389     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(ctab_SSE[4],ctab_SSE[5],ctab_SSE[6],ctab_SSE[7],ctabt_SSE[1],ctabt_SSE[3]); \
 390                                                                         \
 391     GMX_2_MM_TO_M256(ctabt_SSE[0],ctabt_SSE[1],ctab0_SSE);              \
 392     GMX_2_MM_TO_M256(ctabt_SSE[2],ctabt_SSE[3],ctab1_SSE);              \
 393 }
 394
 395 #define load_table_f_v(tab_coul_FDV0, ti_SSE, ti, ctab0_SSE, ctab1_SSE, ctabv_SSE) \
 396 {                                                                       \
 397     __m128 ctab_SSE[8],ctabt_SSE[4],ctabvt_SSE[2];                      \
 398     int    j;                                                           \
 399                                                                         \
 400     /* Bit shifting would be faster, but AVX doesn't support that */    \
 401     _mm256_store_si256((__m256i *)ti,ti_SSE);                           \
 402     for(j=0; j<8; j++)                                                  \
 403     {                                                                   \
 404         ctab_SSE[j] = _mm_load_ps(tab_coul_FDV0+ti[j]*4);               \
 405     }                                                                   \
 406     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(ctab_SSE[0],ctab_SSE[1],ctab_SSE[2],ctab_SSE[3],ctabt_SSE[0],ctabt_SSE[2]); \
 407     GMX_MM_SHUFFLE_4_PS_FIL01_TO_2_PS(ctab_SSE[4],ctab_SSE[5],ctab_SSE[6],ctab_SSE[7],ctabt_SSE[1],ctabt_SSE[3]); \
 408                                                                         \
 409     GMX_2_MM_TO_M256(ctabt_SSE[0],ctabt_SSE[1],ctab0_SSE);              \
 410     GMX_2_MM_TO_M256(ctabt_SSE[2],ctabt_SSE[3],ctab1_SSE);              \
 411                                                                         \
 412     GMX_MM_SHUFFLE_4_PS_FIL2_TO_1_PS(ctab_SSE[0],ctab_SSE[1],ctab_SSE[2],ctab_SSE[3],ctabvt_SSE[0]); \
 413     GMX_MM_SHUFFLE_4_PS_FIL2_TO_1_PS(ctab_SSE[4],ctab_SSE[5],ctab_SSE[6],ctab_SSE[7],ctabvt_SSE[1]); \
 414                                                                         \
 415     GMX_2_MM_TO_M256(ctabvt_SSE[0],ctabvt_SSE[1],ctabv_SSE);            \
 416 }
 417
 418 #endif
 419
 420 #if defined GMX_MM128_HERE && defined GMX_DOUBLE
 421
 422 #define load_table_f(tab_coul_F, ti_SSE, ti, ctab0_SSE, ctab1_SSE)      \
 423 {                                                                       \
 424     int     idx[2];                                                     \
 425     __m128d ctab_SSE[2];                                                \
 426                                                                         \
 427     /* Without SSE4.1 the extract macro needs an immediate: unroll */   \
 428     idx[0] = gmx_mm_extract_epi32(ti_SSE,0);                            \
 429     ctab_SSE[0] = _mm_loadu_pd(tab_coul_F+idx[0]);                      \
 430     idx[1] = gmx_mm_extract_epi32(ti_SSE,1);                            \
 431     ctab_SSE[1] = _mm_loadu_pd(tab_coul_F+idx[1]);                      \
 432                                                                         \
 433     /* Shuffle the force table entries to a convenient order */         \
 434     GMX_MM_TRANSPOSE2_OP_PD(ctab_SSE[0],ctab_SSE[1],ctab0_SSE,ctab1_SSE); \
 435     /* The second force table entry should contain the difference */    \
 436     ctab1_SSE = _mm_sub_pd(ctab1_SSE,ctab0_SSE);                        \
 437 }
 438
 439 #define load_table_f_v(tab_coul_F, tab_coul_V, ti_SSE, ti, ctab0_SSE, ctab1_SSE, ctabv_SSE) \
 440 {                                                                       \
 441     int     idx[2];                                                     \
 442     __m128d ctab_SSE[4];                                                \
 443                                                                         \
 444     /* Without SSE4.1 the extract macro needs an immediate: unroll */   \
 445     idx[0] = gmx_mm_extract_epi32(ti_SSE,0);                            \
 446     ctab_SSE[0] = _mm_loadu_pd(tab_coul_F+idx[0]);                      \
 447     idx[1] = gmx_mm_extract_epi32(ti_SSE,1);                            \
 448     ctab_SSE[1] = _mm_loadu_pd(tab_coul_F+idx[1]);                      \
 449                                                                         \
 450     /* Shuffle the force table entries to a convenient order */         \
 451     GMX_MM_TRANSPOSE2_OP_PD(ctab_SSE[0],ctab_SSE[1],ctab0_SSE,ctab1_SSE); \
 452     /* The second force table entry should contain the difference */    \
 453     ctab1_SSE = _mm_sub_pd(ctab1_SSE,ctab0_SSE);                        \
 454                                                                         \
 455     ctab_SSE[2] = _mm_loadu_pd(tab_coul_V+idx[0]);                      \
 456     ctab_SSE[3] = _mm_loadu_pd(tab_coul_V+idx[1]);                      \
 457                                                                         \
 458     /* Shuffle the energy table entries to a single register */         \
 459     ctabv_SSE = _mm_shuffle_pd(ctab_SSE[2],ctab_SSE[3],_MM_SHUFFLE2(0,0)); \
 460 }
 461
 462 #endif
 463
 464 #if defined GMX_MM256_HERE && defined GMX_DOUBLE
 465
 466 /* Put two 128-bit 2-double registers into one 256-bit 4-ouble register */
 467 #define GMX_2_M128D_TO_M256D(in0,in1,out)                               \
 468 {                                                                       \
 469     out = _mm256_insertf128_pd(_mm256_castpd128_pd256(in0),in1,1);      \
 470 }
 471
 472 #define load_table_f(tab_coul_F, ti_SSE, ti, ctab0_SSE, ctab1_SSE)      \
 473 {                                                                       \
 474     __m128d ctab_SSE[4],tr_SSE[4];                                      \
 475     int     j;                                                          \
 476                                                                         \
 477     _mm_store_si128((__m128i *)ti,ti_SSE);                              \
 478     for(j=0; j<4; j++)                                                  \
 479     {                                                                   \
 480         ctab_SSE[j] = _mm_loadu_pd(tab_coul_F+ti[j]);                   \
 481     }                                                                   \
 482     /* Shuffle the force table entries to a convenient order */         \
 483     GMX_MM_TRANSPOSE2_OP_PD(ctab_SSE[0],ctab_SSE[1],tr_SSE[0],tr_SSE[1]); \
 484     GMX_MM_TRANSPOSE2_OP_PD(ctab_SSE[2],ctab_SSE[3],tr_SSE[2],tr_SSE[3]); \
 485     GMX_2_M128D_TO_M256D(tr_SSE[0],tr_SSE[2],ctab0_SSE);                \
 486     GMX_2_M128D_TO_M256D(tr_SSE[1],tr_SSE[3],ctab1_SSE);                \
 487     /* The second force table entry should contain the difference */    \
 488     ctab1_SSE = _mm256_sub_pd(ctab1_SSE,ctab0_SSE);                     \
 489 }
 490
 491 #define load_table_f_v(tab_coul_F, tab_coul_V, ti_SSE, ti, ctab0_SSE, ctab1_SSE, ctabv_SSE) \
 492 {                                                                       \
 493     __m128d ctab_SSE[8],tr_SSE[4];                                      \
 494     int     j;                                                          \
 495                                                                         \
 496     _mm_store_si128((__m128i *)ti,ti_SSE);                              \
 497     for(j=0; j<4; j++)                                                  \
 498     {                                                                   \
 499         ctab_SSE[j] = _mm_loadu_pd(tab_coul_F+ti[j]);                   \
 500     }                                                                   \
 501     /* Shuffle the force table entries to a convenient order */         \
 502     GMX_MM_TRANSPOSE2_OP_PD(ctab_SSE[0],ctab_SSE[1],tr_SSE[0],tr_SSE[1]); \
 503     GMX_MM_TRANSPOSE2_OP_PD(ctab_SSE[2],ctab_SSE[3],tr_SSE[2],tr_SSE[3]); \
 504     GMX_2_M128D_TO_M256D(tr_SSE[0],tr_SSE[2],ctab0_SSE);                \
 505     GMX_2_M128D_TO_M256D(tr_SSE[1],tr_SSE[3],ctab1_SSE);                \
 506     /* The second force table entry should contain the difference */    \
 507     ctab1_SSE = _mm256_sub_pd(ctab1_SSE,ctab0_SSE);                     \
 508                                                                         \
 509     for(j=0; j<4; j++)                                                  \
 510     {                                                                   \
 511         ctab_SSE[4+j] = _mm_loadu_pd(tab_coul_V+ti[j]);                 \
 512     }                                                                   \
 513     /* Shuffle the energy table entries to a single register */         \
 514     GMX_2_M128D_TO_M256D(_mm_shuffle_pd(ctab_SSE[4],ctab_SSE[5],_MM_SHUFFLE2(0,0)),_mm_shuffle_pd(ctab_SSE[6],ctab_SSE[7],_MM_SHUFFLE2(0,0)),ctabv_SSE); \
 515 }
 516
 517 #endif
 518
 519
 520 /* Add energy register to possibly multiple terms in the energy array.
 521  * This function is the same for SSE/AVX single/double.
 522  */
 523 static inline void add_ener_grp(gmx_mm_pr e_SSE,real *v,const int *offset_jj)
 524 {
 525     int jj;
 526
 527     /* We need to balance the number of store operations with
 528      * the rapidly increases number of combinations of energy groups.
 529      * We add to a temporary buffer for 1 i-group vs 2 j-groups.
 530      */
 531     for(jj=0; jj<(UNROLLJ/2); jj++)
 532     {
 533         gmx_mm_pr v_SSE;
 534
 535         v_SSE = gmx_load_pr(v+offset_jj[jj]+jj*GMX_SIMD_WIDTH_HERE);
 536         gmx_store_pr(v+offset_jj[jj]+jj*GMX_SIMD_WIDTH_HERE,gmx_add_pr(v_SSE,e_SSE));
 537     }
 538 }
 539
 540 #if defined GMX_X86_AVX_256 && GMX_SIMD_WIDTH_HERE == 8
 541 /* As add_ener_grp above, but for two groups of UNROLLJ/2 stored in
 542  * a single SIMD register.
 543  */
 544 static inline void add_ener_grp_halves(gmx_mm_pr e_SSE,
 545                                        real *v0,real *v1,const int *offset_jj)
 546 {
 547     gmx_mm_hpr e_SSE0,e_SSE1;
 548     int jj;
 549
 550     e_SSE0 = _mm256_extractf128_ps(e_SSE,0);
 551     e_SSE1 = _mm256_extractf128_ps(e_SSE,1);
 552
 553     for(jj=0; jj<(UNROLLJ/2); jj++)
 554     {
 555         gmx_mm_hpr v_SSE;
 556
 557         v_SSE = gmx_load_hpr(v0+offset_jj[jj]+jj*GMX_SIMD_WIDTH_HERE/2);
 558         gmx_store_hpr(v0+offset_jj[jj]+jj*GMX_SIMD_WIDTH_HERE/2,gmx_add_hpr(v_SSE,e_SSE0));
 559     }
 560     for(jj=0; jj<(UNROLLJ/2); jj++)
 561     {
 562         gmx_mm_hpr v_SSE;
 563
 564         v_SSE = gmx_load_hpr(v1+offset_jj[jj]+jj*GMX_SIMD_WIDTH_HERE/2);
 565         gmx_store_hpr(v1+offset_jj[jj]+jj*GMX_SIMD_WIDTH_HERE/2,gmx_add_hpr(v_SSE,e_SSE1));
 566     }
 567 }
 568 #endif
 569
 570 #endif /* GMX_X86_SSE2 */
 571
 572 #endif /* _nbnxn_kernel_sse_utils_h_ */