src/gromacs/nbnxm/nbnxm_setup.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2019, by the GROMACS development team, led by
   5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   6  * and including many others, as listed in the AUTHORS file in the
   7  * top-level source directory and at http://www.gromacs.org.
   8  *
   9  * GROMACS is free software; you can redistribute it and/or
  10  * modify it under the terms of the GNU Lesser General Public License
  11  * as published by the Free Software Foundation; either version 2.1
  12  * of the License, or (at your option) any later version.
  13  *
  14  * GROMACS is distributed in the hope that it will be useful,
  15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  17  * Lesser General Public License for more details.
  18  *
  19  * You should have received a copy of the GNU Lesser General Public
  20  * License along with GROMACS; if not, see
  21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  23  *
  24  * If you want to redistribute modifications to GROMACS, please
  25  * consider that scientific software is very special. Version
  26  * control is crucial - bugs must be traceable. We will be happy to
  27  * consider code for inclusion in the official distribution, but
  28  * derived work must not be called official GROMACS. Details are found
  29  * in the README & COPYING files - if they are missing, get the
  30  * official version at http://www.gromacs.org.
  31  *
  32  * To help us fund GROMACS development, we humbly ask that you cite
  33  * the research papers on the package. Check out http://www.gromacs.org.
  34  */
  35 /*! \internal \file
  36  * \brief Common functions for the different NBNXN GPU implementations.
  37  *
  38  * \author Berk Hess <hess@kth.se>
  39  *
  40  * \ingroup module_nbnxm
  41  */
  42
  43 #include "gmxpre.h"
  44
  45 #include "gromacs/domdec/domdec.h"
  46 #include "gromacs/domdec/domdec_struct.h"
  47 #include "gromacs/hardware/hw_info.h"
  48 #include "gromacs/mdlib/gmx_omp_nthreads.h"
  49 #include "gromacs/mdtypes/commrec.h"
  50 #include "gromacs/mdtypes/forcerec.h"
  51 #include "gromacs/mdtypes/inputrec.h"
  52 #include "gromacs/nbnxm/atomdata.h"
  53 #include "gromacs/nbnxm/gpu_data_mgmt.h"
  54 #include "gromacs/nbnxm/nbnxm.h"
  55 #include "gromacs/nbnxm/nbnxm_geometry.h"
  56 #include "gromacs/nbnxm/nbnxm_simd.h"
  57 #include "gromacs/nbnxm/pairlist_tuning.h"
  58 #include "gromacs/simd/simd.h"
  59 #include "gromacs/utility/fatalerror.h"
  60 #include "gromacs/utility/logger.h"
  61
  62 #include "internal.h"
  63
  64 /*! \brief Returns whether CPU SIMD support exists for the given inputrec
  65  *
  66  * If the return value is FALSE and fplog/cr != NULL, prints a fallback
  67  * message to fplog/stderr.
  68  */
  69 static gmx_bool nbnxn_simd_supported(const gmx::MDLogger &mdlog,
  70                                      const t_inputrec    *ir)
  71 {
  72     if (ir->vdwtype == evdwPME && ir->ljpme_combination_rule == eljpmeLB)
  73     {
  74         /* LJ PME with LB combination rule does 7 mesh operations.
  75          * This so slow that we don't compile SIMD non-bonded kernels
  76          * for that. */
  77         GMX_LOG(mdlog.warning).asParagraph().appendText("LJ-PME with Lorentz-Berthelot is not supported with SIMD kernels, falling back to plain C kernels");
  78         return FALSE;
  79     }
  80
  81     return TRUE;
  82 }
  83
  84 /*! \brief Returns the most suitable CPU kernel type and Ewald handling */
  85 static void pick_nbnxn_kernel_cpu(const t_inputrec gmx_unused    *ir,
  86                                   int                            *kernel_type,
  87                                   int                            *ewald_excl,
  88                                   const gmx_hw_info_t gmx_unused &hardwareInfo)
  89 {
  90     *kernel_type = nbnxnk4x4_PlainC;
  91     *ewald_excl  = ewaldexclTable;
  92
  93 #if GMX_SIMD
  94     {
  95 #ifdef GMX_NBNXN_SIMD_4XN
  96         *kernel_type = nbnxnk4xN_SIMD_4xN;
  97 #endif
  98 #ifdef GMX_NBNXN_SIMD_2XNN
  99         *kernel_type = nbnxnk4xN_SIMD_2xNN;
 100 #endif
 101
 102 #if defined GMX_NBNXN_SIMD_2XNN && defined GMX_NBNXN_SIMD_4XN
 103         /* We need to choose if we want 2x(N+N) or 4xN kernels.
 104          * This is based on the SIMD acceleration choice and CPU information
 105          * detected at runtime.
 106          *
 107          * 4xN calculates more (zero) interactions, but has less pair-search
 108          * work and much better kernel instruction scheduling.
 109          *
 110          * Up till now we have only seen that on Intel Sandy/Ivy Bridge,
 111          * which doesn't have FMA, both the analytical and tabulated Ewald
 112          * kernels have similar pair rates for 4x8 and 2x(4+4), so we choose
 113          * 2x(4+4) because it results in significantly fewer pairs.
 114          * For RF, the raw pair rate of the 4x8 kernel is higher than 2x(4+4),
 115          * 10% with HT, 50% without HT. As we currently don't detect the actual
 116          * use of HT, use 4x8 to avoid a potential performance hit.
 117          * On Intel Haswell 4x8 is always faster.
 118          */
 119         *kernel_type = nbnxnk4xN_SIMD_4xN;
 120
 121 #if !GMX_SIMD_HAVE_FMA
 122         if (EEL_PME_EWALD(ir->coulombtype) ||
 123             EVDW_PME(ir->vdwtype))
 124         {
 125             /* We have Ewald kernels without FMA (Intel Sandy/Ivy Bridge).
 126              * There are enough instructions to make 2x(4+4) efficient.
 127              */
 128             *kernel_type = nbnxnk4xN_SIMD_2xNN;
 129         }
 130 #endif
 131         if (hardwareInfo.haveAmdZenCpu)
 132         {
 133             /* One 256-bit FMA per cycle makes 2xNN faster */
 134             *kernel_type = nbnxnk4xN_SIMD_2xNN;
 135         }
 136 #endif  /* GMX_NBNXN_SIMD_2XNN && GMX_NBNXN_SIMD_4XN */
 137
 138
 139         if (getenv("GMX_NBNXN_SIMD_4XN") != nullptr)
 140         {
 141 #ifdef GMX_NBNXN_SIMD_4XN
 142             *kernel_type = nbnxnk4xN_SIMD_4xN;
 143 #else
 144             gmx_fatal(FARGS, "SIMD 4xN kernels requested, but GROMACS has been compiled without support for these kernels");
 145 #endif
 146         }
 147         if (getenv("GMX_NBNXN_SIMD_2XNN") != nullptr)
 148         {
 149 #ifdef GMX_NBNXN_SIMD_2XNN
 150             *kernel_type = nbnxnk4xN_SIMD_2xNN;
 151 #else
 152             gmx_fatal(FARGS, "SIMD 2x(N+N) kernels requested, but GROMACS has been compiled without support for these kernels");
 153 #endif
 154         }
 155
 156         /* Analytical Ewald exclusion correction is only an option in
 157          * the SIMD kernel.
 158          * Since table lookup's don't parallelize with SIMD, analytical
 159          * will probably always be faster for a SIMD width of 8 or more.
 160          * With FMA analytical is sometimes faster for a width if 4 as well.
 161          * In single precision, this is faster on Bulldozer.
 162          */
 163 #if GMX_SIMD_REAL_WIDTH >= 8 || \
 164         (GMX_SIMD_REAL_WIDTH >= 4 && GMX_SIMD_HAVE_FMA && !GMX_DOUBLE)
 165         /* On AMD Zen, tabulated Ewald kernels are faster on all 4 combinations
 166          * of single or double precision and 128 or 256-bit AVX2.
 167          */
 168         if (!hardwareInfo.haveAmdZenCpu)
 169         {
 170             *ewald_excl = ewaldexclAnalytical;
 171         }
 172 #endif
 173         if (getenv("GMX_NBNXN_EWALD_TABLE") != nullptr)
 174         {
 175             *ewald_excl = ewaldexclTable;
 176         }
 177         if (getenv("GMX_NBNXN_EWALD_ANALYTICAL") != nullptr)
 178         {
 179             *ewald_excl = ewaldexclAnalytical;
 180         }
 181
 182     }
 183 #endif // GMX_SIMD
 184 }
 185
 186 const char *lookup_nbnxn_kernel_name(int kernel_type)
 187 {
 188     const char *returnvalue = nullptr;
 189     switch (kernel_type)
 190     {
 191         case nbnxnkNotSet:
 192             returnvalue = "not set";
 193             break;
 194         case nbnxnk4x4_PlainC:
 195             returnvalue = "plain C";
 196             break;
 197         case nbnxnk4xN_SIMD_4xN:
 198         case nbnxnk4xN_SIMD_2xNN:
 199 #if GMX_SIMD
 200             returnvalue = "SIMD";
 201 #else  // GMX_SIMD
 202             returnvalue = "not available";
 203 #endif // GMX_SIMD
 204             break;
 205         case nbnxnk8x8x8_GPU: returnvalue    = "GPU"; break;
 206         case nbnxnk8x8x8_PlainC: returnvalue = "plain C"; break;
 207
 208         case nbnxnkNR:
 209         default:
 210             gmx_fatal(FARGS, "Illegal kernel type selected");
 211     }
 212     return returnvalue;
 213 };
 214
 215 /*! \brief Returns the most suitable kernel type and Ewald handling */
 216 static void pick_nbnxn_kernel(const gmx::MDLogger &mdlog,
 217                               gmx_bool             use_simd_kernels,
 218                               const gmx_hw_info_t &hardwareInfo,
 219                               gmx_bool             bUseGPU,
 220                               EmulateGpuNonbonded  emulateGpu,
 221                               const t_inputrec    *ir,
 222                               int                 *kernel_type,
 223                               int                 *ewald_excl,
 224                               gmx_bool             bDoNonbonded)
 225 {
 226     GMX_RELEASE_ASSERT(kernel_type, "Need a valid kernel_type pointer");
 227
 228     *kernel_type = nbnxnkNotSet;
 229     *ewald_excl  = ewaldexclTable;
 230
 231     if (emulateGpu == EmulateGpuNonbonded::Yes)
 232     {
 233         *kernel_type = nbnxnk8x8x8_PlainC;
 234
 235         if (bDoNonbonded)
 236         {
 237             GMX_LOG(mdlog.warning).asParagraph().appendText("Emulating a GPU run on the CPU (slow)");
 238         }
 239     }
 240     else if (bUseGPU)
 241     {
 242         *kernel_type = nbnxnk8x8x8_GPU;
 243     }
 244
 245     if (*kernel_type == nbnxnkNotSet)
 246     {
 247         if (use_simd_kernels &&
 248             nbnxn_simd_supported(mdlog, ir))
 249         {
 250             pick_nbnxn_kernel_cpu(ir, kernel_type, ewald_excl, hardwareInfo);
 251         }
 252         else
 253         {
 254             *kernel_type = nbnxnk4x4_PlainC;
 255         }
 256     }
 257
 258     if (bDoNonbonded)
 259     {
 260         GMX_LOG(mdlog.info).asParagraph().appendTextFormatted(
 261                 "Using %s %dx%d nonbonded short-range kernels",
 262                 lookup_nbnxn_kernel_name(*kernel_type),
 263                 nbnxn_kernel_to_cluster_i_size(*kernel_type),
 264                 nbnxn_kernel_to_cluster_j_size(*kernel_type));
 265
 266         if (nbnxnk4x4_PlainC == *kernel_type ||
 267             nbnxnk8x8x8_PlainC == *kernel_type)
 268         {
 269             GMX_LOG(mdlog.warning).asParagraph().appendTextFormatted(
 270                     "WARNING: Using the slow %s kernels. This should\n"
 271                     "not happen during routine usage on supported platforms.",
 272                     lookup_nbnxn_kernel_name(*kernel_type));
 273         }
 274     }
 275 }
 276
 277 void init_nb_verlet(const gmx::MDLogger     &mdlog,
 278                     nonbonded_verlet_t     **nb_verlet,
 279                     gmx_bool                 bFEP_NonBonded,
 280                     const t_inputrec        *ir,
 281                     const t_forcerec        *fr,
 282                     const t_commrec         *cr,
 283                     const gmx_hw_info_t     &hardwareInfo,
 284                     const gmx_device_info_t *deviceInfo,
 285                     const gmx_mtop_t        *mtop,
 286                     matrix                   box)
 287 {
 288     nonbonded_verlet_t *nbv;
 289     char               *env;
 290
 291     nbv = new nonbonded_verlet_t();
 292
 293     nbv->emulateGpu = ((getenv("GMX_EMULATE_GPU") != nullptr) ? EmulateGpuNonbonded::Yes : EmulateGpuNonbonded::No);
 294     nbv->bUseGPU    = deviceInfo != nullptr;
 295
 296     GMX_RELEASE_ASSERT(!(nbv->emulateGpu == EmulateGpuNonbonded::Yes && nbv->bUseGPU), "When GPU emulation is active, there cannot be a GPU assignment");
 297
 298     nbv->nbs             = nullptr;
 299     nbv->min_ci_balanced = 0;
 300
 301     nbv->ngrp = (DOMAINDECOMP(cr) ? 2 : 1);
 302     for (int i = 0; i < nbv->ngrp; i++)
 303     {
 304         nbv->grp[i].nbl_lists.nnbl = 0;
 305         nbv->grp[i].kernel_type    = nbnxnkNotSet;
 306
 307         if (i == 0) /* local */
 308         {
 309             pick_nbnxn_kernel(mdlog, fr->use_simd_kernels, hardwareInfo,
 310                               nbv->bUseGPU, nbv->emulateGpu, ir,
 311                               &nbv->grp[i].kernel_type,
 312                               &nbv->grp[i].ewald_excl,
 313                               fr->bNonbonded);
 314         }
 315         else /* non-local */
 316         {
 317             /* Use the same kernel for local and non-local interactions */
 318             nbv->grp[i].kernel_type = nbv->grp[0].kernel_type;
 319             nbv->grp[i].ewald_excl  = nbv->grp[0].ewald_excl;
 320         }
 321     }
 322
 323     nbv->listParams = std::make_unique<NbnxnListParameters>(ir->rlist);
 324     setupDynamicPairlistPruning(mdlog, ir, mtop, box, nbv->grp[0].kernel_type, fr->ic,
 325                                 nbv->listParams.get());
 326
 327     nbv->nbs = std::make_unique<nbnxn_search>(DOMAINDECOMP(cr) ? &cr->dd->nc : nullptr,
 328                                               DOMAINDECOMP(cr) ? domdec_zones(cr->dd) : nullptr,
 329                                               bFEP_NonBonded,
 330                                               gmx_omp_nthreads_get(emntPairsearch));
 331
 332     for (int i = 0; i < nbv->ngrp; i++)
 333     {
 334         nbnxn_init_pairlist_set(&nbv->grp[i].nbl_lists,
 335                                 nbnxn_kernel_pairlist_simple(nbv->grp[i].kernel_type),
 336                                 /* 8x8x8 "non-simple" lists are ATM always combined */
 337                                 !nbnxn_kernel_pairlist_simple(nbv->grp[i].kernel_type));
 338     }
 339
 340     int      enbnxninitcombrule;
 341     if (fr->ic->vdwtype == evdwCUT &&
 342         (fr->ic->vdw_modifier == eintmodNONE ||
 343          fr->ic->vdw_modifier == eintmodPOTSHIFT) &&
 344         getenv("GMX_NO_LJ_COMB_RULE") == nullptr)
 345     {
 346         /* Plain LJ cut-off: we can optimize with combination rules */
 347         enbnxninitcombrule = enbnxninitcombruleDETECT;
 348     }
 349     else if (fr->ic->vdwtype == evdwPME)
 350     {
 351         /* LJ-PME: we need to use a combination rule for the grid */
 352         if (fr->ljpme_combination_rule == eljpmeGEOM)
 353         {
 354             enbnxninitcombrule = enbnxninitcombruleGEOM;
 355         }
 356         else
 357         {
 358             enbnxninitcombrule = enbnxninitcombruleLB;
 359         }
 360     }
 361     else
 362     {
 363         /* We use a full combination matrix: no rule required */
 364         enbnxninitcombrule = enbnxninitcombruleNONE;
 365     }
 366
 367     nbv->nbat = new nbnxn_atomdata_t(nbv->bUseGPU ? gmx::PinningPolicy::PinnedIfSupported : gmx::PinningPolicy::CannotBePinned);
 368     int mimimumNumEnergyGroupNonbonded = ir->opts.ngener;
 369     if (ir->opts.ngener - ir->nwall == 1)
 370     {
 371         /* We have only one non-wall energy group, we do not need energy group
 372          * support in the non-bondeds kernels, since all non-bonded energy
 373          * contributions go to the first element of the energy group matrix.
 374          */
 375         mimimumNumEnergyGroupNonbonded = 1;
 376     }
 377     bool bSimpleList = nbnxn_kernel_pairlist_simple(nbv->grp[0].kernel_type);
 378     nbnxn_atomdata_init(mdlog,
 379                         nbv->nbat,
 380                         nbv->grp[0].kernel_type,
 381                         enbnxninitcombrule,
 382                         fr->ntype, fr->nbfp,
 383                         mimimumNumEnergyGroupNonbonded,
 384                         bSimpleList ? gmx_omp_nthreads_get(emntNonbonded) : 1);
 385
 386     if (nbv->bUseGPU)
 387     {
 388         /* init the NxN GPU data; the last argument tells whether we'll have
 389          * both local and non-local NB calculation on GPU */
 390         nbnxn_gpu_init(&nbv->gpu_nbv,
 391                        deviceInfo,
 392                        fr->ic,
 393                        nbv->listParams.get(),
 394                        nbv->nbat,
 395                        cr->nodeid,
 396                        (nbv->ngrp > 1));
 397
 398         if ((env = getenv("GMX_NB_MIN_CI")) != nullptr)
 399         {
 400             char *end;
 401
 402             nbv->min_ci_balanced = strtol(env, &end, 10);
 403             if (!end || (*end != 0) || nbv->min_ci_balanced < 0)
 404             {
 405                 gmx_fatal(FARGS, "Invalid value passed in GMX_NB_MIN_CI=%s, non-negative integer required", env);
 406             }
 407
 408             if (debug)
 409             {
 410                 fprintf(debug, "Neighbor-list balancing parameter: %d (passed as env. var.)\n",
 411                         nbv->min_ci_balanced);
 412             }
 413         }
 414         else
 415         {
 416             nbv->min_ci_balanced = nbnxn_gpu_min_ci_balanced(nbv->gpu_nbv);
 417             if (debug)
 418             {
 419                 fprintf(debug, "Neighbor-list balancing parameter: %d (auto-adjusted to the number of GPU multi-processors)\n",
 420                         nbv->min_ci_balanced);
 421             }
 422         }
 423
 424     }
 425
 426     *nb_verlet = nbv;
 427 }