src/gromacs/domdec/domdec_internal.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2014,2015,2016,2017,2018 by the GROMACS development team.
   5  * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
   6  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   7  * and including many others, as listed in the AUTHORS file in the
   8  * top-level source directory and at http://www.gromacs.org.
   9  *
  10  * GROMACS is free software; you can redistribute it and/or
  11  * modify it under the terms of the GNU Lesser General Public License
  12  * as published by the Free Software Foundation; either version 2.1
  13  * of the License, or (at your option) any later version.
  14  *
  15  * GROMACS is distributed in the hope that it will be useful,
  16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  18  * Lesser General Public License for more details.
  19  *
  20  * You should have received a copy of the GNU Lesser General Public
  21  * License along with GROMACS; if not, see
  22  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  23  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  24  *
  25  * If you want to redistribute modifications to GROMACS, please
  26  * consider that scientific software is very special. Version
  27  * control is crucial - bugs must be traceable. We will be happy to
  28  * consider code for inclusion in the official distribution, but
  29  * derived work must not be called official GROMACS. Details are found
  30  * in the README & COPYING files - if they are missing, get the
  31  * official version at http://www.gromacs.org.
  32  *
  33  * To help us fund GROMACS development, we humbly ask that you cite
  34  * the research papers on the package. Check out http://www.gromacs.org.
  35  */
  36 /*! \internal \file
  37  *
  38  * \brief Declares implementation functions and types for the domain
  39  * decomposition module.
  40  *
  41  * \author Berk Hess <hess@kth.se>
  42  * \ingroup module_domdec
  43  */
  44 #ifndef GMX_DOMDEC_DOMDEC_INTERNAL_H
  45 #define GMX_DOMDEC_DOMDEC_INTERNAL_H
  46
  47 #include "config.h"
  48
  49 #include "gromacs/domdec/domdec.h"
  50 #include "gromacs/domdec/domdec_struct.h"
  51 #include "gromacs/mdlib/updategroupscog.h"
  52 #include "gromacs/timing/cyclecounter.h"
  53 #include "gromacs/topology/block.h"
  54
  55 struct t_commrec;
  56
  57 /*! \cond INTERNAL */
  58
  59 #define DD_NLOAD_MAX 9
  60
  61 struct BalanceRegion;
  62
  63 namespace gmx
  64 {
  65 enum class DdRankOrder : int;
  66 }
  67 // namespace
  68
  69
  70 //! Indices to communicate in a dimension
  71 struct gmx_domdec_ind_t
  72 {
  73     //! @{
  74     /*! \brief The numbers of charge groups to send and receive for each
  75      * cell that requires communication, the last entry contains the total
  76      * number of atoms that needs to be communicated.
  77      */
  78     int nsend[DD_MAXIZONE + 2] = {};
  79     int nrecv[DD_MAXIZONE + 2] = {};
  80     //! @}
  81     //! The charge groups to send
  82     std::vector<int> index;
  83     //! @{
  84     /* The atom range for non-in-place communication */
  85     int cell2at0[DD_MAXIZONE] = {};
  86     int cell2at1[DD_MAXIZONE] = {};
  87     //! @}
  88 };
  89
  90 //! Things relating to index communication
  91 struct gmx_domdec_comm_dim_t
  92 {
  93     /* Returns the number of grid pulses (the number of domains in the halo along this dimension) */
  94     int numPulses() const { return ind.size(); }
  95
  96     /**< For dlb, for use with edlbAUTO          */
  97     int np_dlb = 0;
  98     /**< The indices to communicate, size np     */
  99     std::vector<gmx_domdec_ind_t> ind;
 100     /**< Can we receive data in place?            */
 101     bool receiveInPlace = false;
 102 };
 103
 104 /*! \brief Load balancing data along a dim used on the master rank of that dim */
 105 struct RowMaster
 106 {
 107     struct Bounds
 108     {
 109         /**< State var.: max lower bound., incl. neighbors */
 110         real cellFracLowerMax = 0;
 111         /**< State var.: min upper bound., incl. neighbors */
 112         real cellFracUpperMin = 0;
 113         /**< Temp. var.: lower limit for cell boundary */
 114         real boundMin = 0;
 115         /**< Temp. var.: upper limit for cell boundary */
 116         real boundMax = 0;
 117     };
 118
 119     /**< Temp. var.: is this cell size at the limit */
 120     std::vector<bool> isCellMin;
 121     /**< State var.: cell boundaries, box relative */
 122     std::vector<real> cellFrac;
 123     /**< Temp. var.: old cell size */
 124     std::vector<real> oldCellFrac;
 125     /**< Cell bounds */
 126     std::vector<Bounds> bounds;
 127     /**< State var.: is DLB limited in this row */
 128     bool dlbIsLimited = false;
 129     /**< Temp. var.  */
 130     std::vector<real> buf_ncd;
 131 };
 132
 133 /*! \brief Struct for managing cell sizes with DLB along a dimension */
 134 struct DDCellsizesWithDlb
 135 {
 136     /**< Cell row root struct, only available on the first rank in a row */
 137     std::unique_ptr<RowMaster> rowMaster;
 138     /**< The cell sizes, in fractions, along a row, not available on the first rank in a row */
 139     std::vector<real> fracRow;
 140     /**< The lower corner, in fractions, in triclinic space */
 141     real fracLower = 0;
 142     /**< The upper corner, in fractions, in triclinic space */
 143     real fracUpper = 0;
 144     /**< The maximum lower corner among all our neighbors */
 145     real fracLowerMax = 0;
 146     /**< The minimum upper corner among all our neighbors */
 147     real fracUpperMin = 0;
 148 };
 149
 150 /*! \brief Struct for compute load commuication
 151  *
 152  * Here floats are accurate enough, since these variables
 153  * only influence the load balancing, not the actual MD results.
 154  */
 155 typedef struct domdec_load
 156 {
 157     /**< The number of load recordings */
 158     int nload = 0;
 159     /**< Scan of the sum of load over dimensions */
 160     float* load = nullptr;
 161     /**< The sum of the load over the ranks up to our current dimension */
 162     float sum = 0;
 163     /**< The maximum over the ranks contributing to \p sum */
 164     float max = 0;
 165     /**< Like \p sum, but takes the maximum when the load balancing is limited */
 166     float sum_m = 0;
 167     /**< Minimum cell volume, relative to the box */
 168     float cvol_min = 0;
 169     /**< The PP time during which PME can overlap */
 170     float mdf = 0;
 171     /**< The PME-only rank load */
 172     float pme = 0;
 173     /**< Bit flags that tell if DLB was limited, per dimension */
 174     int flags = 0;
 175 } domdec_load_t;
 176
 177 /*! \brief Data needed to sort an atom to the desired location in the local state */
 178 typedef struct gmx_cgsort
 179 {
 180     /**< Neighborsearch grid cell index */
 181     int nsc = 0;
 182     /**< Global atom/charge group index */
 183     int ind_gl = 0;
 184     /**< Local atom/charge group index */
 185     int ind = 0;
 186 } gmx_cgsort_t;
 187
 188 /*! \brief Temporary buffers for sorting atoms */
 189 typedef struct gmx_domdec_sort
 190 {
 191     /**< Sorted array of indices */
 192     std::vector<gmx_cgsort_t> sorted;
 193     /**< Array of stationary atom/charge group indices */
 194     std::vector<gmx_cgsort_t> stationary;
 195     /**< Array of moved atom/charge group indices */
 196     std::vector<gmx_cgsort_t> moved;
 197     /**< Integer buffer for sorting */
 198     std::vector<int> intBuffer;
 199     /**< Int64 buffer for sorting */
 200     std::vector<int64_t> int64Buffer;
 201 } gmx_domdec_sort_t;
 202
 203 /*! \brief Manages atom ranges and order for the local state atom vectors */
 204 class DDAtomRanges
 205 {
 206 public:
 207     /*! \brief The local state atom order
 208      *
 209      * This enum determines the order of the atoms in the local state.
 210      * ddnatHOME and ddnatZONE should be first and second,
 211      * the others can be ordered as wanted.
 212      */
 213     enum class Type : int
 214     {
 215         Home,        /**< The home atoms */
 216         Zones,       /**< All zones in the eighth shell */
 217         Vsites,      /**< Atoms communicated for virtual sites */
 218         Constraints, /**< Atoms communicated for constraints */
 219         Number       /**< Not a count, only present for convenience */
 220     };
 221
 222     /*! \brief Returns the start atom index for range \p rangeType */
 223     int start(Type rangeType) const
 224     {
 225         if (rangeType == Type::Home)
 226         {
 227             return 0;
 228         }
 229         else
 230         {
 231             return end_[static_cast<int>(rangeType) - 1];
 232         }
 233     }
 234
 235     /*! \brief Returns the end atom index for range \p rangeType */
 236     int end(Type rangeType) const { return end_[static_cast<int>(rangeType)]; }
 237
 238     /*! \brief Returns the number of home atoms */
 239     int numHomeAtoms() const { return end_[static_cast<int>(Type::Home)]; }
 240
 241     /*! \brief Returns the total number of atoms */
 242     int numAtomsTotal() const { return end_[static_cast<int>(Type::Number) - 1]; }
 243
 244     /*! \brief Sets the end index of range \p rangeType to \p end
 245      *
 246      * This should be called either with Type::Home or with a type
 247      * that is larger than that passed in the previous call to setEnd.
 248      * A release assertion for these conditions are present.
 249      */
 250     void setEnd(Type rangeType, int end)
 251     {
 252         GMX_RELEASE_ASSERT(rangeType == Type::Home || rangeType > lastTypeSet_,
 253                            "Can only set either home or a larger type than the last one");
 254
 255         for (int i = static_cast<int>(rangeType); i < static_cast<int>(Type::Number); i++)
 256         {
 257             end_[i] = end;
 258         }
 259
 260         lastTypeSet_ = rangeType;
 261     }
 262
 263 private:
 264     /*! \brief The list of end atom indices */
 265     std::array<int, static_cast<int>(Type::Number)> end_;
 266     Type                                            lastTypeSet_ = Type::Number;
 267 };
 268
 269 /*! \brief Enum of dynamic load balancing states
 270  *
 271  * Allowed DLB states and transitions
 272  * - intialization at startup:
 273  *                             -> offUser ("-dlb no")
 274  *                             -> onUser  ("-dlb yes")
 275  *                             -> offCanTurnOn ("-dlb auto")
 276  *
 277  * - in automatic mode (i.e. initial state offCanTurnOn):
 278  *   offCanTurnOn         -> onCanTurnOff
 279  *   offCanTurnOn         -> offForever
 280  *   offCanTurnOn         -> offTemporarilyLocked
 281  *   offTemporarilyLocked -> offCanTurnOn
 282  *   onCanTurnOff         -> offCanTurnOn
 283  */
 284 enum class DlbState
 285 {
 286     offUser,    /**< DLB is permanently off per user request */
 287     offForever, /**< DLB is off due to a runtime condition (not supported or causes performance loss) and will never be turned on */
 288     offCanTurnOn,         /**< DLB is off and will turn on on imbalance */
 289     offTemporarilyLocked, /**< DLB is off and temporarily can't turn on */
 290     onCanTurnOff,         /**< DLB is on and can turn off when slow */
 291     onUser,               /**< DLB is permanently on per user request */
 292     Count                 /**< The number of DLB states */
 293 };
 294
 295 /*! \brief The PME domain decomposition for one dimension */
 296 typedef struct gmx_ddpme
 297 {
 298     /**< The dimension */
 299     int dim = 0;
 300     /**< Tells if DD and PME dims match */
 301     gmx_bool dim_match = false;
 302     /**< The number of PME ranks/domains in this dimension */
 303     int nslab = 0;
 304     /**< Cell sizes for determining the PME comm. with SLB */
 305     real* slb_dim_f = nullptr;
 306     /**< The minimum pp node location, size nslab */
 307     int* pp_min = nullptr;
 308     /**< The maximum pp node location, size nslab */
 309     int* pp_max = nullptr;
 310     /**< The maximum shift for coordinate redistribution in PME */
 311     int maxshift = 0;
 312 } gmx_ddpme_t;
 313
 314 struct gmx_ddzone_t
 315 {
 316     /**< The minimum bottom of this zone                        */
 317     real min0 = 0;
 318     /**< The maximum top of this zone                           */
 319     real max1 = 0;
 320     /**< The minimum top of this zone                           */
 321     real min1 = 0;
 322     /**< The maximum bottom communicaton height for this zone   */
 323     real mch0 = 0;
 324     /**< The maximum top communicaton height for this zone      */
 325     real mch1 = 0;
 326     /**< The bottom value of the first cell in this zone        */
 327     real p1_0 = 0;
 328     /**< The top value of the first cell in this zone           */
 329     real p1_1 = 0;
 330     /**< Bool disguised as a real, 1 when the above data has been set. 0 otherwise */
 331     real dataSet = 0;
 332 };
 333
 334 /*! \brief The number of reals in gmx_ddzone_t */
 335 constexpr int c_ddzoneNumReals = 8;
 336
 337 template<typename T>
 338 class DDBufferAccess;
 339
 340 /*! \brief Temporary storage container that minimizes (re)allocation and clearing
 341  *
 342  * This is only the storage, actual access happens through DDBufferAccess.
 343  * All methods check if the buffer is (not) in use.
 344  */
 345 template<typename T>
 346 class DDBuffer
 347 {
 348 private:
 349     /*! \brief Returns a buffer of size \p numElements, the elements are undefined */
 350     gmx::ArrayRef<T> resize(size_t numElements)
 351     {
 352         GMX_ASSERT(isInUse_, "Should only operate on acquired buffers");
 353
 354         if (numElements > buffer_.size())
 355         {
 356             buffer_.resize(numElements);
 357         }
 358
 359         return gmx::arrayRefFromArray(buffer_.data(), numElements);
 360     }
 361
 362     /*! \brief Acquire the buffer for use with size set to \p numElements, the elements are undefined */
 363     gmx::ArrayRef<T> acquire(size_t numElements)
 364     {
 365         GMX_RELEASE_ASSERT(!isInUse_, "Should only request free buffers");
 366         isInUse_ = true;
 367
 368         return resize(numElements);
 369     }
 370
 371     /*! \brief Releases the buffer, buffer_ should not be used after this */
 372     void release()
 373     {
 374         GMX_RELEASE_ASSERT(isInUse_, "Should only release buffers in use");
 375         isInUse_ = false;
 376     }
 377
 378     std::vector<T> buffer_;          /**< The actual memory buffer */
 379     bool           isInUse_ = false; /**< Flag that tells whether the buffer is in use */
 380
 381     friend class DDBufferAccess<T>;
 382 };
 383
 384 /*! \brief Class that manages access to a temporary memory buffer */
 385 template<typename T>
 386 class DDBufferAccess
 387 {
 388 public:
 389     /*! \brief Constructor, returns a buffer of size \p numElements, element values are undefined
 390      *
 391      * \note The actual memory buffer \p ddBuffer can not be used to
 392      *       create other DDBufferAccess objects until the one created
 393      *       here is destroyed.
 394      */
 395     DDBufferAccess(DDBuffer<T>& ddBuffer, size_t numElements) : ddBuffer_(ddBuffer)
 396     {
 397         buffer = ddBuffer_.acquire(numElements);
 398     }
 399
 400     ~DDBufferAccess() { ddBuffer_.release(); }
 401
 402     /*! \brief Resizes the buffer to \p numElements, new elements are undefined
 403      *
 404      * \note The buffer arrayref is updated after this call.
 405      */
 406     void resize(size_t numElements) { buffer = ddBuffer_.resize(numElements); }
 407
 408 private:
 409     DDBuffer<T>& ddBuffer_; /**< Reference to the storage class */
 410 public:
 411     gmx::ArrayRef<T> buffer; /**< The access to the memory buffer */
 412 };
 413
 414 /*! \brief Temporary buffer for setting up communiation over one pulse and all zones in the halo */
 415 struct dd_comm_setup_work_t
 416 {
 417     /**< The local atom group indices to send */
 418     std::vector<int> localAtomGroupBuffer;
 419     /**< Buffer for collecting the global atom group indices to send */
 420     std::vector<int> atomGroupBuffer;
 421     /**< Buffer for collecting the atom group positions to send */
 422     std::vector<gmx::RVec> positionBuffer;
 423     /**< The number of atoms contained in the atom groups to send */
 424     int nat = 0;
 425     /**< The number of atom groups to send for the last zone */
 426     int nsend_zone = 0;
 427 };
 428
 429 /*! \brief Information about the simulated system */
 430 struct DDSystemInfo
 431 {
 432     //! True when update groups are used
 433     bool useUpdateGroups = false;
 434     //! Update atom grouping for each molecule type
 435     gmx::ArrayRef<const gmx::RangePartitioning> updateGroupingsPerMoleculeType;
 436     //! The maximum radius over all update groups
 437     real maxUpdateGroupRadius;
 438
 439     //! Are molecules always whole, i.e. not broken by PBC?
 440     bool moleculesAreAlwaysWhole = false;
 441     //! Are there inter-domain bonded interactions?
 442     bool haveInterDomainBondeds = false;
 443     //! Are there inter-domain multi-body interactions?
 444     bool haveInterDomainMultiBodyBondeds = false;
 445
 446     //! Cut-off for multi-body interactions
 447     real minCutoffForMultiBody = 0;
 448     //! Cut-off for non-bonded/2-body interactions
 449     real cutoff = 0;
 450     //! The lower limit for the DD cell size
 451     real cellsizeLimit = 0;
 452
 453     //! Can atoms connected by constraints be assigned to different domains?
 454     bool mayHaveSplitConstraints = false;
 455     //! Can atoms connected by settles be assigned to different domains?
 456     bool mayHaveSplitSettles = false;
 457     //! Estimated communication range needed for constraints
 458     real constraintCommunicationRange = 0;
 459
 460     //! Whether to only communicate atoms beyond the non-bonded cut-off when they are involved in bonded interactions with non-local atoms
 461     bool filterBondedCommunication = false;
 462     //! Whether to increase the multi-body cut-off beyond the minimum required
 463     bool increaseMultiBodyCutoff = false;
 464 };
 465
 466 /*! \brief Settings that affect the behavior of the domain decomposition
 467  *
 468  * These settings depend on options chosen by the user, set by enviroment
 469  * variables, as well as hardware support. The initial DLB state also
 470  * depends on the integrator.
 471  *
 472  * Note: Settings that depend on the simulated system are in DDSystemInfo.
 473  */
 474 struct DDSettings
 475 {
 476     //! Use MPI_Sendrecv communication instead of non-blocking calls
 477     bool useSendRecv2 = false;
 478
 479     /* Information for managing the dynamic load balancing */
 480     //! Maximum DLB scaling per load balancing step in percent
 481     int dlb_scale_lim = 0;
 482     //! Flop counter (0=no,1=yes,2=with (eFlop-1)*5% noise
 483     int eFlop = 0;
 484
 485     //! Whether to order the DD dimensions from z to x
 486     bool useDDOrderZYX = false;
 487
 488     //! Whether to use MPI Cartesian reordering of communicators, when supported (almost never)
 489     bool useCartesianReorder = true;
 490
 491     //! Whether we should record the load
 492     bool recordLoad = false;
 493
 494     /* Debugging */
 495     //! Step interval for dumping the local+non-local atoms to pdb
 496     int nstDDDump = 0;
 497     //! Step interval for duming the DD grid to pdb
 498     int nstDDDumpGrid = 0;
 499     //! DD debug print level: 0, 1, 2
 500     int DD_debug = 0;
 501
 502     //! The DLB state at the start of the run
 503     DlbState initialDlbState = DlbState::offCanTurnOn;
 504 };
 505
 506 /*! \brief Information on how the DD ranks are set up */
 507 struct DDRankSetup
 508 {
 509     /**< The rank ordering */
 510     gmx::DdRankOrder rankOrder;
 511
 512     /**< The number of particle-particle (non PME-only) ranks */
 513     int numPPRanks = 0;
 514     /**< The DD PP grid */
 515     ivec numPPCells = { 0, 0, 0 };
 516
 517     /* PME and Cartesian communicator stuff */
 518     bool usePmeOnlyRanks = false;
 519     /**< The number of decomposition dimensions for PME, 0: no PME */
 520     int npmedecompdim = 0;
 521     /**< The number of ranks doing PME (PP/PME or only PME) */
 522     int numRanksDoingPme = 0;
 523     /**< The number of PME ranks/domains along x */
 524     int npmenodes_x = 0;
 525     /**< The number of PME ranks/domains along y */
 526     int npmenodes_y = 0;
 527     /**< The 1D or 2D PME domain decomposition setup */
 528     gmx_ddpme_t ddpme[2];
 529 };
 530
 531 /*! \brief Information on Cartesian MPI setup of the DD ranks */
 532 struct CartesianRankSetup
 533 {
 534     /**< Use Cartesian communication between PP and PME ranks */
 535     bool bCartesianPP_PME = false;
 536     /**< Cartesian grid for combinted PP+PME ranks */
 537     ivec ntot = {};
 538     /**< The number of dimensions for the PME setup that are Cartesian */
 539     int cartpmedim = 0;
 540     /**< The Cartesian index to sim rank conversion, used with bCartesianPP_PME */
 541     std::vector<int> ddindex2simnodeid;
 542
 543     /* The DD particle-particle nodes only */
 544     /**< Use a Cartesian communicator for PP */
 545     bool bCartesianPP = false;
 546     /**< The Cartesian index to DD rank conversion, used with bCartesianPP */
 547     std::vector<int> ddindex2ddnodeid;
 548 };
 549
 550 /*! \brief Struct for domain decomposition communication
 551  *
 552  * This struct contains most information about domain decomposition
 553  * communication setup, some communication buffers, some statistics
 554  * and also the setup for the communication between particle-particle
 555  * and PME only ranks.
 556  *
 557  * All arrays are indexed with 0 to dd->ndim (not Cartesian indexing),
 558  * unless stated otherwise.
 559  */
 560 struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
 561 {
 562     /**< Constant parameters that control DD behavior */
 563     DDSettings ddSettings;
 564
 565     /**< Information on how the DD ranks are set up */
 566     DDRankSetup ddRankSetup;
 567     /**< Information on the Cartesian part of the DD rank setup */
 568     CartesianRankSetup cartesianRankSetup;
 569
 570     /* Charge group / atom sorting */
 571     /**< Data structure for cg/atom sorting */
 572     std::unique_ptr<gmx_domdec_sort_t> sort;
 573
 574     //! Centers of mass of local update groups
 575     std::unique_ptr<gmx::UpdateGroupsCog> updateGroupsCog;
 576
 577     /* Data for the optional filtering of communication of atoms for bonded interactions */
 578     /**< Links between atoms through bonded interactions */
 579     t_blocka* bondedLinks = nullptr;
 580
 581     /* The DLB state, possible values are defined above */
 582     DlbState dlbState;
 583     /* With dlbState=DlbState::offCanTurnOn, should we check if to DLB on at the next DD? */
 584     gmx_bool bCheckWhetherToTurnDlbOn = false;
 585     /* The first DD count since we are running without DLB */
 586     int ddPartioningCountFirstDlbOff = 0;
 587
 588     /* Cell sizes for static load balancing, first index cartesian */
 589     real** slb_frac = nullptr;
 590
 591     /**< Information about the simulated system */
 592     DDSystemInfo systemInfo;
 593
 594     /* The width of the communicated boundaries */
 595     /**< Cut-off for multi-body interactions, also 2-body bonded when \p cutoff_mody > \p cutoff */
 596     real cutoff_mbody = 0;
 597     /**< The minimum guaranteed cell-size, Cartesian indexing */
 598     gmx::RVec cellsize_min = { 0, 0, 0 };
 599     /**< The minimum guaranteed cell-size with dlb=auto */
 600     gmx::RVec cellsize_min_dlb = { 0, 0, 0 };
 601     /**< The lower limit for the DD cell size with DLB */
 602     real cellsize_limit = 0;
 603     /**< Effectively no NB cut-off limit with DLB for systems without PBC? */
 604     bool bVacDLBNoLimit = false;
 605
 606     /** With PME load balancing we set limits on DLB */
 607     bool bPMELoadBalDLBLimits = false;
 608     /** DLB needs to take into account that we want to allow this maximum
 609      *  cut-off (for PME load balancing), this could limit cell boundaries.
 610      */
 611     real PMELoadBal_max_cutoff = 0;
 612
 613     /**< box lower corner, required with dim's without pbc and -gcom */
 614     gmx::RVec box0 = { 0, 0, 0 };
 615     /**< box size, required with dim's without pbc and -gcom */
 616     gmx::RVec box_size = { 0, 0, 0 };
 617
 618     /**< The DD cell lower corner, in triclinic space */
 619     gmx::RVec cell_x0 = { 0, 0, 0 };
 620     /**< The DD cell upper corner, in triclinic space */
 621     gmx::RVec cell_x1 = { 0, 0, 0 };
 622
 623     /**< The old \p cell_x0, to check cg displacements */
 624     gmx::RVec old_cell_x0 = { 0, 0, 0 };
 625     /**< The old \p cell_x1, to check cg displacements */
 626     gmx::RVec old_cell_x1 = { 0, 0, 0 };
 627
 628     /** The communication setup and charge group boundaries for the zones */
 629     gmx_domdec_zones_t zones;
 630
 631     /* The zone limits for DD dimensions 1 and 2 (not 0), determined from
 632      * cell boundaries of neighboring cells for staggered grids when using
 633      * dynamic load balancing.
 634      */
 635     /**< Zone limits for dim 1 with staggered grids */
 636     std::array<gmx_ddzone_t, 2> zone_d1;
 637     /**< Zone limits for dim 2 with staggered grids */
 638     gmx_ddzone_t zone_d2[2][2];
 639
 640     /** The coordinate/force communication setup and indices */
 641     std::array<gmx_domdec_comm_dim_t, DIM> cd;
 642     /** Restricts the maximum number of cells to communicate with in one dimension
 643      *
 644      * Dynamic load balancing is not permitted to change sizes if it
 645      * would violate this restriction. */
 646     int maxpulse = 0;
 647
 648     /** Which cg distribution is stored on the master node,
 649      *  stored as DD partitioning call count.
 650      */
 651     int64_t master_cg_ddp_count = 0;
 652
 653     /** The number of cg's received from the direct neighbors */
 654     std::array<int, DD_MAXZONE> zone_ncg1 = { 0 };
 655
 656     /** The atom ranges in the local state */
 657     DDAtomRanges atomRanges;
 658
 659     /** Array for signalling if atoms have moved to another domain */
 660     std::vector<int> movedBuffer;
 661
 662     /** Communication int buffer for general use */
 663     DDBuffer<int> intBuffer;
 664
 665     /** Communication rvec buffer for general use */
 666     DDBuffer<gmx::RVec> rvecBuffer;
 667
 668     /* Temporary storage for thread parallel communication setup */
 669     /**< Thread-local work data */
 670     std::vector<dd_comm_setup_work_t> dth;
 671
 672     /* Communication buffer only used with multiple grid pulses */
 673     /**< Another rvec comm. buffer */
 674     DDBuffer<gmx::RVec> rvecBuffer2;
 675
 676     /* Communication buffers for local redistribution */
 677     /**< Charge group flag comm. buffers */
 678     std::array<std::vector<int>, DIM * 2> cggl_flag;
 679     /**< Charge group center comm. buffers */
 680     std::array<std::vector<gmx::RVec>, DIM * 2> cgcm_state;
 681
 682     /* Cell sizes for dynamic load balancing */
 683     std::vector<DDCellsizesWithDlb> cellsizesWithDlb;
 684
 685     /* Stuff for load communication */
 686     /**< The recorded load data */
 687     domdec_load_t* load = nullptr;
 688     /**< The number of MPI ranks sharing the GPU our rank is using */
 689     int nrank_gpu_shared = 0;
 690 #if GMX_MPI
 691     /**< The MPI load communicator */
 692     MPI_Comm* mpi_comm_load = nullptr;
 693     /**< The MPI load communicator for ranks sharing a GPU */
 694     MPI_Comm mpi_comm_gpu_shared;
 695 #endif
 696
 697     /**< Struct for timing the force load balancing region */
 698     BalanceRegion* balanceRegion = nullptr;
 699
 700     /* Cycle counters over nstlist steps */
 701     /**< Total cycles counted */
 702     std::array<float, ddCyclNr> cycl = { 0 };
 703     /**< The number of cycle recordings */
 704     std::array<int, ddCyclNr> cycl_n = { 0 };
 705     /**< The maximum cycle count */
 706     std::array<float, ddCyclNr> cycl_max = { 0 };
 707     /**< Total flops counted */
 708     double flop = 0.0;
 709     /**< The number of flop recordings */
 710     int flop_n = 0;
 711     /** How many times did we have load measurements */
 712     int n_load_have = 0;
 713     /** How many times have we collected the load measurements */
 714     int n_load_collect = 0;
 715
 716     /* Cycle count history for DLB checks */
 717     /**< The averaged cycles per step over the last nstlist step before turning on DLB */
 718     float cyclesPerStepBeforeDLB = 0;
 719     /**< The running average of the cycles per step during DLB */
 720     float cyclesPerStepDlbExpAverage = 0;
 721     /**< Have we turned off DLB (after turning DLB on)? */
 722     bool haveTurnedOffDlb = false;
 723     /**< The DD step at which we last measured that DLB off was faster than DLB on, 0 if there was no such step */
 724     int64_t dlbSlowerPartitioningCount = 0;
 725
 726     /* Statistics for atoms */
 727     /**< The atoms per range, summed over the steps */
 728     double sum_nat[static_cast<int>(DDAtomRanges::Type::Number)] = {};
 729
 730     /* Statistics for calls and times */
 731     /**< The number of partioning calls */
 732     int ndecomp = 0;
 733     /**< The number of load recordings */
 734     int nload = 0;
 735     /**< Total MD step time */
 736     double load_step = 0.0;
 737     /**< Total PP force time */
 738     double load_sum = 0.0;
 739     /**< Max \p load_sum over the ranks */
 740     double load_max = 0.0;
 741     /**< Was load balancing limited, per DD dim */
 742     gmx::IVec load_lim = { 0, 0, 0 };
 743     /**< Total time on PP done during PME overlap time */
 744     double load_mdf = 0.0;
 745     /**< Total time on our PME-only rank */
 746     double load_pme = 0.0;
 747
 748     /** The last partition step */
 749     int64_t partition_step = INT_MIN;
 750 };
 751
 752 /*! \brief DD zone permutation
 753  *
 754  * Zone permutation from the Cartesian x-major/z-minor order to an order
 755  * that leads to consecutive charge groups for neighbor searching.
 756  * TODO: It should be possible to remove this now that the group scheme is removed
 757  */
 758 static const int zone_perm[3][4] = { { 0, 0, 0, 0 }, { 1, 0, 0, 0 }, { 3, 0, 1, 2 } };
 759
 760 /*! \brief DD zone reordering to Cartesian order
 761  *
 762  * Index to reorder the zone such that the end up in Cartesian order
 763  * with dimension index 0 major and dimension index 2 minor.
 764  */
 765 static const int zone_reorder_cartesian[DD_MAXZONE] = { 0, 1, 3, 2, 5, 4, 6, 7 };
 766
 767 /* dd_zo and dd_zp3 is set up such that i zones with non-zero
 768  * components see only j zones with that component 0.
 769  */
 770
 771 /*! \brief Returns the DD cut-off distance for multi-body interactions */
 772 real dd_cutoff_multibody(const gmx_domdec_t* dd);
 773
 774 /*! \brief Returns the DD cut-off distance for two-body interactions */
 775 real dd_cutoff_twobody(const gmx_domdec_t* dd);
 776
 777 /*! \brief Returns the domain index given the number of domains and the domain coordinates
 778  *
 779  * This order is required to minimize the coordinate communication in PME
 780  * which uses decomposition in the x direction.
 781  */
 782 static inline int dd_index(const ivec numDomains, const ivec domainCoordinates)
 783 {
 784     return ((domainCoordinates[XX] * numDomains[YY] + domainCoordinates[YY]) * numDomains[ZZ])
 785            + domainCoordinates[ZZ];
 786 };
 787
 788 /*! Returns the size of the buffer to hold fractional cell boundaries for DD dimension index dimIndex */
 789 static inline int ddCellFractionBufferSize(const gmx_domdec_t* dd, int dimIndex)
 790 {
 791     return dd->numCells[dd->dim[dimIndex]] + 1 + dimIndex * 2 + 1 + dimIndex;
 792 }
 793
 794 /*! \brief Maximum number of ranks for using send/recv for state scattering and gathering
 795  *
 796  * Use separate MPI send and receive commands
 797  * when #ranks <= c_maxNumRanksUseSendRecvForScatterAndGather
 798  * This saves memory (and some copying for small #ranks).
 799  * For high parallelization scatter and gather calls are used.
 800  */
 801 static constexpr int c_maxNumRanksUseSendRecvForScatterAndGather = 4;
 802
 803 /*! \brief Make DD cells larger by this factor than the limit to avoid rounding issues */
 804 static constexpr double DD_CELL_MARGIN = 1.0001;
 805
 806 /*! \brief Factor for checking DD cell size limitation during DLB, should be in between 1 and DD_CELL_MARGIN */
 807 static constexpr double DD_CELL_MARGIN2 = 1.00005;
 808
 809 /*! \brief With pressure scaling, keep cell sizes 2% above the limit to allow for some scaling */
 810 static constexpr double DD_PRES_SCALE_MARGIN = 1.02;
 811
 812 /*! \endcond */
 813
 814 #endif