src/mdlib/nbnxn_kernels/nbnxn_kernel_simd_4xn.c

   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  * check out http://www.gromacs.org for more information.
   7  * Copyright (c) 2012,2013, by the GROMACS development team, led by
   8  * David van der Spoel, Berk Hess, Erik Lindahl, and including many
   9  * others, as listed in the AUTHORS file in the top-level source
  10  * directory and at http://www.gromacs.org.
  11  *
  12  * GROMACS is free software; you can redistribute it and/or
  13  * modify it under the terms of the GNU Lesser General Public License
  14  * as published by the Free Software Foundation; either version 2.1
  15  * of the License, or (at your option) any later version.
  16  *
  17  * GROMACS is distributed in the hope that it will be useful,
  18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  20  * Lesser General Public License for more details.
  21  *
  22  * You should have received a copy of the GNU Lesser General Public
  23  * License along with GROMACS; if not, see
  24  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  25  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  26  *
  27  * If you want to redistribute modifications to GROMACS, please
  28  * consider that scientific software is very special. Version
  29  * control is crucial - bugs must be traceable. We will be happy to
  30  * consider code for inclusion in the official distribution, but
  31  * derived work must not be called official GROMACS. Details are found
  32  * in the README & COPYING files - if they are missing, get the
  33  * official version at http://www.gromacs.org.
  34  *
  35  * To help us fund GROMACS development, we humbly ask that you cite
  36  * the research papers on the package. Check out http://www.gromacs.org.
  37  */
  38 #ifdef HAVE_CONFIG_H
  39 #include <config.h>
  40 #endif
  41
  42 #include <math.h>
  43
  44 #include "typedefs.h"
  45 #include "vec.h"
  46 #include "smalloc.h"
  47 #include "force.h"
  48 #include "gmx_omp_nthreads.h"
  49 #include "../nbnxn_consts.h"
  50 #include "nbnxn_kernel_common.h"
  51
  52 #ifdef GMX_NBNXN_SIMD_4XN
  53
  54 #include "nbnxn_kernel_simd_4xn.h"
  55
  56 /* Include all flavors of the SSE or AVX 4xN kernel loops */
  57
  58 #if GMX_NBNXN_SIMD_BITWIDTH == 128
  59 #define GMX_MM128_HERE
  60 #else
  61 #if GMX_NBNXN_SIMD_BITWIDTH == 256
  62 #define GMX_MM256_HERE
  63 #else
  64 #error "unsupported GMX_NBNXN_SIMD_BITWIDTH"
  65 #endif
  66 #endif
  67
  68 /* Analytical reaction-field kernels */
  69 #define CALC_COUL_RF
  70
  71 #include "nbnxn_kernel_simd_4xn_includes.h"
  72
  73 #undef CALC_COUL_RF
  74
  75 /* Tabulated exclusion interaction electrostatics kernels */
  76 #define CALC_COUL_TAB
  77
  78 /* Single cut-off: rcoulomb = rvdw */
  79 #include "nbnxn_kernel_simd_4xn_includes.h"
  80
  81 /* Twin cut-off: rcoulomb >= rvdw */
  82 #define VDW_CUTOFF_CHECK
  83 #include "nbnxn_kernel_simd_4xn_includes.h"
  84 #undef VDW_CUTOFF_CHECK
  85
  86 #undef CALC_COUL_TAB
  87
  88 /* Analytical Ewald exclusion interaction electrostatics kernels */
  89 #define CALC_COUL_EWALD
  90
  91 /* Single cut-off: rcoulomb = rvdw */
  92 #include "nbnxn_kernel_simd_4xn_includes.h"
  93
  94 /* Twin cut-off: rcoulomb >= rvdw */
  95 #define VDW_CUTOFF_CHECK
  96 #include "nbnxn_kernel_simd_4xn_includes.h"
  97 #undef VDW_CUTOFF_CHECK
  98
  99 #undef CALC_COUL_EWALD
 100
 101
 102 typedef void (*p_nbk_func_ener)(const nbnxn_pairlist_t     *nbl,
 103                                 const nbnxn_atomdata_t     *nbat,
 104                                 const interaction_const_t  *ic,
 105                                 rvec                       *shift_vec,
 106                                 real                       *f,
 107                                 real                       *fshift,
 108                                 real                       *Vvdw,
 109                                 real                       *Vc);
 110
 111 typedef void (*p_nbk_func_noener)(const nbnxn_pairlist_t     *nbl,
 112                                   const nbnxn_atomdata_t     *nbat,
 113                                   const interaction_const_t  *ic,
 114                                   rvec                       *shift_vec,
 115                                   real                       *f,
 116                                   real                       *fshift);
 117
 118 enum {
 119     coultRF, coultTAB, coultTAB_TWIN, coultEWALD, coultEWALD_TWIN, coultNR
 120 };
 121
 122 #define NBK_FN(elec, ljcomb) nbnxn_kernel_simd_4xn_ ## elec ## _comb_ ## ljcomb ## _ener
 123 static p_nbk_func_ener p_nbk_ener[coultNR][ljcrNR] =
 124 { { NBK_FN(rf, geom), NBK_FN(rf, lb), NBK_FN(rf, none) },
 125   { NBK_FN(tab, geom), NBK_FN(tab, lb), NBK_FN(tab, none) },
 126   { NBK_FN(tab_twin, geom), NBK_FN(tab_twin, lb), NBK_FN(tab_twin, none) },
 127   { NBK_FN(ewald, geom), NBK_FN(ewald, lb), NBK_FN(ewald, none) },
 128   { NBK_FN(ewald_twin, geom), NBK_FN(ewald_twin, lb), NBK_FN(ewald_twin, none) } };
 129 #undef NBK_FN
 130
 131 #define NBK_FN(elec, ljcomb) nbnxn_kernel_simd_4xn_ ## elec ## _comb_ ## ljcomb ## _energrp
 132 static p_nbk_func_ener p_nbk_energrp[coultNR][ljcrNR] =
 133 { { NBK_FN(rf, geom), NBK_FN(rf, lb), NBK_FN(rf, none) },
 134   { NBK_FN(tab, geom), NBK_FN(tab, lb), NBK_FN(tab, none) },
 135   { NBK_FN(tab_twin, geom), NBK_FN(tab_twin, lb), NBK_FN(tab_twin, none) },
 136   { NBK_FN(ewald, geom), NBK_FN(ewald, lb), NBK_FN(ewald, none) },
 137   { NBK_FN(ewald_twin, geom), NBK_FN(ewald_twin, lb), NBK_FN(ewald_twin, none) } };
 138 #undef NBK_FN
 139
 140 #define NBK_FN(elec, ljcomb) nbnxn_kernel_simd_4xn_ ## elec ## _comb_ ## ljcomb ## _noener
 141 static p_nbk_func_noener p_nbk_noener[coultNR][ljcrNR] =
 142 { { NBK_FN(rf, geom), NBK_FN(rf, lb), NBK_FN(rf, none) },
 143   { NBK_FN(tab, geom), NBK_FN(tab, lb), NBK_FN(tab, none) },
 144   { NBK_FN(tab_twin, geom), NBK_FN(tab_twin, lb), NBK_FN(tab_twin, none) },
 145   { NBK_FN(ewald, geom), NBK_FN(ewald, lb), NBK_FN(ewald, none) },
 146   { NBK_FN(ewald_twin, geom), NBK_FN(ewald_twin, lb), NBK_FN(ewald_twin, none) } };
 147 #undef NBK_FN
 148
 149
 150 static void reduce_group_energies(int ng, int ng_2log,
 151                                   const real *VSvdw, const real *VSc,
 152                                   real *Vvdw, real *Vc)
 153 {
 154     const int simd_width   = GMX_SIMD_WIDTH_HERE;
 155     const int unrollj_half = GMX_SIMD_WIDTH_HERE/2;
 156     int       ng_p2, i, j, j0, j1, c, s;
 157
 158     ng_p2 = (1<<ng_2log);
 159
 160     /* The size of the x86 SIMD energy group buffer array is:
 161      * ng*ng*ng_p2*unrollj_half*simd_width
 162      */
 163     for (i = 0; i < ng; i++)
 164     {
 165         for (j = 0; j < ng; j++)
 166         {
 167             Vvdw[i*ng+j] = 0;
 168             Vc[i*ng+j]   = 0;
 169         }
 170
 171         for (j1 = 0; j1 < ng; j1++)
 172         {
 173             for (j0 = 0; j0 < ng; j0++)
 174             {
 175                 c = ((i*ng + j1)*ng_p2 + j0)*unrollj_half*simd_width;
 176                 for (s = 0; s < unrollj_half; s++)
 177                 {
 178                     Vvdw[i*ng+j0] += VSvdw[c+0];
 179                     Vvdw[i*ng+j1] += VSvdw[c+1];
 180                     Vc  [i*ng+j0] += VSc  [c+0];
 181                     Vc  [i*ng+j1] += VSc  [c+1];
 182                     c             += simd_width + 2;
 183                 }
 184             }
 185         }
 186     }
 187 }
 188
 189 #endif /* GMX_NBNXN_SIMD_4XN */
 190
 191 void
 192 nbnxn_kernel_simd_4xn(nbnxn_pairlist_set_t       *nbl_list,
 193                       const nbnxn_atomdata_t     *nbat,
 194                       const interaction_const_t  *ic,
 195                       int                         ewald_excl,
 196                       rvec                       *shift_vec,
 197                       int                         force_flags,
 198                       int                         clearF,
 199                       real                       *fshift,
 200                       real                       *Vc,
 201                       real                       *Vvdw)
 202 #ifdef GMX_NBNXN_SIMD_4XN
 203 {
 204     int                nnbl;
 205     nbnxn_pairlist_t **nbl;
 206     int                coult;
 207     int                nb;
 208
 209     nnbl = nbl_list->nnbl;
 210     nbl  = nbl_list->nbl;
 211
 212     if (EEL_RF(ic->eeltype) || ic->eeltype == eelCUT)
 213     {
 214         coult = coultRF;
 215     }
 216     else
 217     {
 218         if (ewald_excl == ewaldexclTable)
 219         {
 220             if (ic->rcoulomb == ic->rvdw)
 221             {
 222                 coult = coultTAB;
 223             }
 224             else
 225             {
 226                 coult = coultTAB_TWIN;
 227             }
 228         }
 229         else
 230         {
 231             if (ic->rcoulomb == ic->rvdw)
 232             {
 233                 coult = coultEWALD;
 234             }
 235             else
 236             {
 237                 coult = coultEWALD_TWIN;
 238             }
 239         }
 240     }
 241
 242 #pragma omp parallel for schedule(static) num_threads(gmx_omp_nthreads_get(emntNonbonded))
 243     for (nb = 0; nb < nnbl; nb++)
 244     {
 245         nbnxn_atomdata_output_t *out;
 246         real                    *fshift_p;
 247
 248         out = &nbat->out[nb];
 249
 250         if (clearF == enbvClearFYes)
 251         {
 252             clear_f(nbat, nb, out->f);
 253         }
 254
 255         if ((force_flags & GMX_FORCE_VIRIAL) && nnbl == 1)
 256         {
 257             fshift_p = fshift;
 258         }
 259         else
 260         {
 261             fshift_p = out->fshift;
 262
 263             if (clearF == enbvClearFYes)
 264             {
 265                 clear_fshift(fshift_p);
 266             }
 267         }
 268
 269         /* With Ewald type electrostatics we the forces for excluded atom pairs
 270          * should not contribute to the virial sum. The exclusion forces
 271          * are not calculate in the energy kernels, but are in _noener.
 272          */
 273         if (!((force_flags & GMX_FORCE_ENERGY) ||
 274               (EEL_FULL(ic->eeltype) && (force_flags & GMX_FORCE_VIRIAL))))
 275         {
 276             /* Don't calculate energies */
 277             p_nbk_noener[coult][nbat->comb_rule](nbl[nb], nbat,
 278                                                  ic,
 279                                                  shift_vec,
 280                                                  out->f,
 281                                                  fshift_p);
 282         }
 283         else if (out->nV == 1 || !(force_flags & GMX_FORCE_ENERGY))
 284         {
 285             /* No energy groups */
 286             out->Vvdw[0] = 0;
 287             out->Vc[0]   = 0;
 288
 289             p_nbk_ener[coult][nbat->comb_rule](nbl[nb], nbat,
 290                                                ic,
 291                                                shift_vec,
 292                                                out->f,
 293                                                fshift_p,
 294                                                out->Vvdw,
 295                                                out->Vc);
 296         }
 297         else
 298         {
 299             /* Calculate energy group contributions */
 300             int i;
 301
 302             for (i = 0; i < out->nVS; i++)
 303             {
 304                 out->VSvdw[i] = 0;
 305             }
 306             for (i = 0; i < out->nVS; i++)
 307             {
 308                 out->VSc[i] = 0;
 309             }
 310
 311             p_nbk_energrp[coult][nbat->comb_rule](nbl[nb], nbat,
 312                                                   ic,
 313                                                   shift_vec,
 314                                                   out->f,
 315                                                   fshift_p,
 316                                                   out->VSvdw,
 317                                                   out->VSc);
 318
 319             reduce_group_energies(nbat->nenergrp, nbat->neg_2log,
 320                                   out->VSvdw, out->VSc,
 321                                   out->Vvdw, out->Vc);
 322         }
 323     }
 324
 325     if (force_flags & GMX_FORCE_ENERGY)
 326     {
 327         reduce_energies_over_lists(nbat, nnbl, Vvdw, Vc);
 328     }
 329 }
 330 #else
 331 {
 332     gmx_incons("nbnxn_kernel_simd_4xn called while GROMACS was configured without 4xN SIMD kernels enabled");
 333 }
 334 #endif