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.h"
  58 #include "gromacs/nbnxm/pairlist_tuning.h"
  59 #include "gromacs/nbnxm/pairlistset.h"
  60 #include "gromacs/simd/simd.h"
  61 #include "gromacs/utility/fatalerror.h"
  62 #include "gromacs/utility/logger.h"
  63
  64 #include "grid.h"
  65 #include "internal.h"
  66
  67 namespace Nbnxm
  68 {
  69
  70 /*! \brief Returns whether CPU SIMD support exists for the given inputrec
  71  *
  72  * If the return value is FALSE and fplog/cr != NULL, prints a fallback
  73  * message to fplog/stderr.
  74  */
  75 static gmx_bool nbnxn_simd_supported(const gmx::MDLogger &mdlog,
  76                                      const t_inputrec    *ir)
  77 {
  78     if (ir->vdwtype == evdwPME && ir->ljpme_combination_rule == eljpmeLB)
  79     {
  80         /* LJ PME with LB combination rule does 7 mesh operations.
  81          * This so slow that we don't compile SIMD non-bonded kernels
  82          * for that. */
  83         GMX_LOG(mdlog.warning).asParagraph().appendText("LJ-PME with Lorentz-Berthelot is not supported with SIMD kernels, falling back to plain C kernels");
  84         return FALSE;
  85     }
  86
  87     return TRUE;
  88 }
  89
  90 /*! \brief Returns the most suitable CPU kernel type and Ewald handling */
  91 static KernelSetup
  92 pick_nbnxn_kernel_cpu(const t_inputrec gmx_unused    *ir,
  93                       const gmx_hw_info_t gmx_unused &hardwareInfo)
  94 {
  95     KernelSetup kernelSetup;
  96
  97     if (!GMX_SIMD)
  98     {
  99         kernelSetup.kernelType         = KernelType::Cpu4x4_PlainC;
 100         kernelSetup.ewaldExclusionType = EwaldExclusionType::Table;
 101     }
 102     else
 103     {
 104 #ifdef GMX_NBNXN_SIMD_4XN
 105         kernelSetup.kernelType = KernelType::Cpu4xN_Simd_4xN;
 106 #endif
 107 #ifdef GMX_NBNXN_SIMD_2XNN
 108         kernelSetup.kernelType = KernelType::Cpu4xN_Simd_2xNN;
 109 #endif
 110
 111 #if defined GMX_NBNXN_SIMD_2XNN && defined GMX_NBNXN_SIMD_4XN
 112         /* We need to choose if we want 2x(N+N) or 4xN kernels.
 113          * This is based on the SIMD acceleration choice and CPU information
 114          * detected at runtime.
 115          *
 116          * 4xN calculates more (zero) interactions, but has less pair-search
 117          * work and much better kernel instruction scheduling.
 118          *
 119          * Up till now we have only seen that on Intel Sandy/Ivy Bridge,
 120          * which doesn't have FMA, both the analytical and tabulated Ewald
 121          * kernels have similar pair rates for 4x8 and 2x(4+4), so we choose
 122          * 2x(4+4) because it results in significantly fewer pairs.
 123          * For RF, the raw pair rate of the 4x8 kernel is higher than 2x(4+4),
 124          * 10% with HT, 50% without HT. As we currently don't detect the actual
 125          * use of HT, use 4x8 to avoid a potential performance hit.
 126          * On Intel Haswell 4x8 is always faster.
 127          */
 128         kernelSetup.kernelType = KernelType::Cpu4xN_Simd_4xN;
 129
 130         if (!GMX_SIMD_HAVE_FMA && (EEL_PME_EWALD(ir->coulombtype) ||
 131                                    EVDW_PME(ir->vdwtype)))
 132         {
 133             /* We have Ewald kernels without FMA (Intel Sandy/Ivy Bridge).
 134              * There are enough instructions to make 2x(4+4) efficient.
 135              */
 136             kernelSetup.kernelType = KernelType::Cpu4xN_Simd_2xNN;
 137         }
 138
 139         if (hardwareInfo.haveAmdZenCpu)
 140         {
 141             /* One 256-bit FMA per cycle makes 2xNN faster */
 142             kernelSetup.kernelType = KernelType::Cpu4xN_Simd_2xNN;
 143         }
 144 #endif      /* GMX_NBNXN_SIMD_2XNN && GMX_NBNXN_SIMD_4XN */
 145
 146
 147         if (getenv("GMX_NBNXN_SIMD_4XN") != nullptr)
 148         {
 149 #ifdef GMX_NBNXN_SIMD_4XN
 150             kernelSetup.kernelType = KernelType::Cpu4xN_Simd_4xN;
 151 #else
 152             gmx_fatal(FARGS, "SIMD 4xN kernels requested, but GROMACS has been compiled without support for these kernels");
 153 #endif
 154         }
 155         if (getenv("GMX_NBNXN_SIMD_2XNN") != nullptr)
 156         {
 157 #ifdef GMX_NBNXN_SIMD_2XNN
 158             kernelSetup.kernelType = KernelType::Cpu4xN_Simd_2xNN;
 159 #else
 160             gmx_fatal(FARGS, "SIMD 2x(N+N) kernels requested, but GROMACS has been compiled without support for these kernels");
 161 #endif
 162         }
 163
 164         /* Analytical Ewald exclusion correction is only an option in
 165          * the SIMD kernel.
 166          * Since table lookup's don't parallelize with SIMD, analytical
 167          * will probably always be faster for a SIMD width of 8 or more.
 168          * With FMA analytical is sometimes faster for a width if 4 as well.
 169          * In single precision, this is faster on Bulldozer.
 170          * On AMD Zen, tabulated Ewald kernels are faster on all 4 combinations
 171          * of single or double precision and 128 or 256-bit AVX2.
 172          */
 173         if (
 174 #if GMX_SIMD
 175             (GMX_SIMD_REAL_WIDTH >= 8 ||
 176              (GMX_SIMD_REAL_WIDTH >= 4 && GMX_SIMD_HAVE_FMA && !GMX_DOUBLE)) &&
 177 #endif
 178             !hardwareInfo.haveAmdZenCpu)
 179         {
 180             kernelSetup.ewaldExclusionType = EwaldExclusionType::Analytical;
 181         }
 182         else
 183         {
 184             kernelSetup.ewaldExclusionType = EwaldExclusionType::Table;
 185         }
 186         if (getenv("GMX_NBNXN_EWALD_TABLE") != nullptr)
 187         {
 188             kernelSetup.ewaldExclusionType = EwaldExclusionType::Table;
 189         }
 190         if (getenv("GMX_NBNXN_EWALD_ANALYTICAL") != nullptr)
 191         {
 192             kernelSetup.ewaldExclusionType = EwaldExclusionType::Analytical;
 193         }
 194
 195     }
 196
 197     return kernelSetup;
 198 }
 199
 200 const char *lookup_kernel_name(const KernelType kernelType)
 201 {
 202     const char *returnvalue = nullptr;
 203     switch (kernelType)
 204     {
 205         case KernelType::NotSet:
 206             returnvalue = "not set";
 207             break;
 208         case KernelType::Cpu4x4_PlainC:
 209             returnvalue = "plain C";
 210             break;
 211         case KernelType::Cpu4xN_Simd_4xN:
 212         case KernelType::Cpu4xN_Simd_2xNN:
 213 #if GMX_SIMD
 214             returnvalue = "SIMD";
 215 #else  // GMX_SIMD
 216             returnvalue = "not available";
 217 #endif // GMX_SIMD
 218             break;
 219         case KernelType::Gpu8x8x8: returnvalue        = "GPU"; break;
 220         case KernelType::Cpu8x8x8_PlainC: returnvalue = "plain C"; break;
 221
 222         default:
 223             gmx_fatal(FARGS, "Illegal kernel type selected");
 224     }
 225     return returnvalue;
 226 };
 227
 228 /*! \brief Returns the most suitable kernel type and Ewald handling */
 229 static KernelSetup
 230 pick_nbnxn_kernel(const gmx::MDLogger     &mdlog,
 231                   gmx_bool                 use_simd_kernels,
 232                   const gmx_hw_info_t     &hardwareInfo,
 233                   const NonbondedResource &nonbondedResource,
 234                   const t_inputrec        *ir,
 235                   gmx_bool                 bDoNonbonded)
 236 {
 237     KernelSetup kernelSetup;
 238
 239     if (nonbondedResource == NonbondedResource::EmulateGpu)
 240     {
 241         kernelSetup.kernelType         = KernelType::Cpu8x8x8_PlainC;
 242         kernelSetup.ewaldExclusionType = EwaldExclusionType::DecidedByGpuModule;
 243
 244         if (bDoNonbonded)
 245         {
 246             GMX_LOG(mdlog.warning).asParagraph().appendText("Emulating a GPU run on the CPU (slow)");
 247         }
 248     }
 249     else if (nonbondedResource == NonbondedResource::Gpu)
 250     {
 251         kernelSetup.kernelType         = KernelType::Gpu8x8x8;
 252         kernelSetup.ewaldExclusionType = EwaldExclusionType::DecidedByGpuModule;
 253     }
 254     else
 255     {
 256         if (use_simd_kernels &&
 257             nbnxn_simd_supported(mdlog, ir))
 258         {
 259             kernelSetup = pick_nbnxn_kernel_cpu(ir, hardwareInfo);
 260         }
 261         else
 262         {
 263             kernelSetup.kernelType         = KernelType::Cpu4x4_PlainC;
 264             kernelSetup.ewaldExclusionType = EwaldExclusionType::Analytical;
 265         }
 266     }
 267
 268     if (bDoNonbonded)
 269     {
 270         GMX_LOG(mdlog.info).asParagraph().appendTextFormatted(
 271                 "Using %s %dx%d nonbonded short-range kernels",
 272                 lookup_kernel_name(kernelSetup.kernelType),
 273                 IClusterSizePerKernelType[kernelSetup.kernelType],
 274                 JClusterSizePerKernelType[kernelSetup.kernelType]);
 275
 276         if (KernelType::Cpu4x4_PlainC == kernelSetup.kernelType ||
 277             KernelType::Cpu8x8x8_PlainC == kernelSetup.kernelType)
 278         {
 279             GMX_LOG(mdlog.warning).asParagraph().appendTextFormatted(
 280                     "WARNING: Using the slow %s kernels. This should\n"
 281                     "not happen during routine usage on supported platforms.",
 282                     lookup_kernel_name(kernelSetup.kernelType));
 283         }
 284     }
 285
 286     GMX_RELEASE_ASSERT(kernelSetup.kernelType != KernelType::NotSet &&
 287                        kernelSetup.ewaldExclusionType != EwaldExclusionType::NotSet,
 288                        "All kernel setup parameters should be set here");
 289
 290     return kernelSetup;
 291 }
 292
 293 } // namespace Nbnxm
 294
 295 void nonbonded_verlet_t::initPairlistSets(const bool haveMultipleDomains)
 296 {
 297     pairlistSets_.emplace_back(*listParams);
 298     if (haveMultipleDomains)
 299     {
 300         pairlistSets_.emplace_back(*listParams);
 301     }
 302 }
 303
 304 namespace Nbnxm
 305 {
 306
 307 void init_nb_verlet(const gmx::MDLogger     &mdlog,
 308                     nonbonded_verlet_t     **nb_verlet,
 309                     gmx_bool                 bFEP_NonBonded,
 310                     const t_inputrec        *ir,
 311                     const t_forcerec        *fr,
 312                     const t_commrec         *cr,
 313                     const gmx_hw_info_t     &hardwareInfo,
 314                     const gmx_device_info_t *deviceInfo,
 315                     const gmx_mtop_t        *mtop,
 316                     matrix                   box)
 317 {
 318     nonbonded_verlet_t *nbv        = new nonbonded_verlet_t();
 319
 320     const bool          emulateGpu = (getenv("GMX_EMULATE_GPU") != nullptr);
 321     const bool          useGpu     = deviceInfo != nullptr;
 322
 323     GMX_RELEASE_ASSERT(!(emulateGpu && useGpu), "When GPU emulation is active, there cannot be a GPU assignment");
 324
 325     NonbondedResource nonbondedResource;
 326     if (useGpu)
 327     {
 328         nonbondedResource = NonbondedResource::Gpu;
 329     }
 330     else if (emulateGpu)
 331     {
 332         nonbondedResource = NonbondedResource::EmulateGpu;
 333     }
 334     else
 335     {
 336         nonbondedResource = NonbondedResource::Cpu;
 337     }
 338
 339     nbv->nbs             = nullptr;
 340
 341     nbv->setKernelSetup(pick_nbnxn_kernel(mdlog, fr->use_simd_kernels, hardwareInfo,
 342                                           nonbondedResource, ir,
 343                                           fr->bNonbonded));
 344
 345     const bool haveMultipleDomains = (DOMAINDECOMP(cr) && cr->dd->nnodes > 1);
 346
 347     nbv->listParams = std::make_unique<NbnxnListParameters>(nbv->kernelSetup().kernelType,
 348                                                             ir->rlist);
 349     nbv->initPairlistSets(haveMultipleDomains);
 350
 351     nbv->min_ci_balanced = 0;
 352
 353     setupDynamicPairlistPruning(mdlog, ir, mtop, box, fr->ic,
 354                                 nbv->listParams.get());
 355
 356     nbv->nbs = std::make_unique<nbnxn_search>(ir->ePBC,
 357                                               DOMAINDECOMP(cr) ? &cr->dd->nc : nullptr,
 358                                               DOMAINDECOMP(cr) ? domdec_zones(cr->dd) : nullptr,
 359                                               bFEP_NonBonded,
 360                                               gmx_omp_nthreads_get(emntPairsearch));
 361
 362     int      enbnxninitcombrule;
 363     if (fr->ic->vdwtype == evdwCUT &&
 364         (fr->ic->vdw_modifier == eintmodNONE ||
 365          fr->ic->vdw_modifier == eintmodPOTSHIFT) &&
 366         getenv("GMX_NO_LJ_COMB_RULE") == nullptr)
 367     {
 368         /* Plain LJ cut-off: we can optimize with combination rules */
 369         enbnxninitcombrule = enbnxninitcombruleDETECT;
 370     }
 371     else if (fr->ic->vdwtype == evdwPME)
 372     {
 373         /* LJ-PME: we need to use a combination rule for the grid */
 374         if (fr->ljpme_combination_rule == eljpmeGEOM)
 375         {
 376             enbnxninitcombrule = enbnxninitcombruleGEOM;
 377         }
 378         else
 379         {
 380             enbnxninitcombrule = enbnxninitcombruleLB;
 381         }
 382     }
 383     else
 384     {
 385         /* We use a full combination matrix: no rule required */
 386         enbnxninitcombrule = enbnxninitcombruleNONE;
 387     }
 388
 389     nbv->nbat = new nbnxn_atomdata_t(useGpu ? gmx::PinningPolicy::PinnedIfSupported : gmx::PinningPolicy::CannotBePinned);
 390     int mimimumNumEnergyGroupNonbonded = ir->opts.ngener;
 391     if (ir->opts.ngener - ir->nwall == 1)
 392     {
 393         /* We have only one non-wall energy group, we do not need energy group
 394          * support in the non-bondeds kernels, since all non-bonded energy
 395          * contributions go to the first element of the energy group matrix.
 396          */
 397         mimimumNumEnergyGroupNonbonded = 1;
 398     }
 399     nbnxn_atomdata_init(mdlog,
 400                         nbv->nbat,
 401                         nbv->kernelSetup().kernelType,
 402                         enbnxninitcombrule,
 403                         fr->ntype, fr->nbfp,
 404                         mimimumNumEnergyGroupNonbonded,
 405                         nbv->pairlistIsSimple() ? gmx_omp_nthreads_get(emntNonbonded) : 1);
 406
 407     if (useGpu)
 408     {
 409         /* init the NxN GPU data; the last argument tells whether we'll have
 410          * both local and non-local NB calculation on GPU */
 411         gpu_init(&nbv->gpu_nbv,
 412                  deviceInfo,
 413                  fr->ic,
 414                  nbv->listParams.get(),
 415                  nbv->nbat,
 416                  cr->nodeid,
 417                  haveMultipleDomains);
 418
 419         if (const char *env = getenv("GMX_NB_MIN_CI"))
 420         {
 421             char *end;
 422
 423             nbv->min_ci_balanced = strtol(env, &end, 10);
 424             if (!end || (*end != 0) || nbv->min_ci_balanced < 0)
 425             {
 426                 gmx_fatal(FARGS, "Invalid value passed in GMX_NB_MIN_CI=%s, non-negative integer required", env);
 427             }
 428
 429             if (debug)
 430             {
 431                 fprintf(debug, "Neighbor-list balancing parameter: %d (passed as env. var.)\n",
 432                         nbv->min_ci_balanced);
 433             }
 434         }
 435         else
 436         {
 437             nbv->min_ci_balanced = gpu_min_ci_balanced(nbv->gpu_nbv);
 438             if (debug)
 439             {
 440                 fprintf(debug, "Neighbor-list balancing parameter: %d (auto-adjusted to the number of GPU multi-processors)\n",
 441                         nbv->min_ci_balanced);
 442             }
 443         }
 444
 445     }
 446
 447     *nb_verlet = nbv;
 448 }
 449
 450 } // namespace Nbnxm