2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2014,2015,2016,2017,2018,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.
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.
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.
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.
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.
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.
37 * \brief Declares implementation functions and types for the domain
38 * decomposition module.
40 * \author Berk Hess <hess@kth.se>
41 * \ingroup module_domdec
43 #ifndef GMX_DOMDEC_DOMDEC_INTERNAL_H
44 #define GMX_DOMDEC_DOMDEC_INTERNAL_H
48 #include "gromacs/domdec/domdec.h"
49 #include "gromacs/domdec/domdec_struct.h"
50 #include "gromacs/mdlib/updategroupscog.h"
51 #include "gromacs/timing/cyclecounter.h"
52 #include "gromacs/topology/block.h"
58 #define DD_NLOAD_MAX 9
62 //! Indices to communicate in a dimension
63 struct gmx_domdec_ind_t
66 /*! \brief The numbers of charge groups to send and receive for each
67 * cell that requires communication, the last entry contains the total
68 * number of atoms that needs to be communicated.
70 int nsend[DD_MAXIZONE + 2] = {};
71 int nrecv[DD_MAXIZONE + 2] = {};
73 //! The charge groups to send
74 std::vector<int> index;
76 /* The atom range for non-in-place communication */
77 int cell2at0[DD_MAXIZONE] = {};
78 int cell2at1[DD_MAXIZONE] = {};
82 //! Things relating to index communication
83 struct gmx_domdec_comm_dim_t
85 /* Returns the number of grid pulses (the number of domains in the halo along this dimension) */
86 int numPulses() const { return ind.size(); }
88 /**< For dlb, for use with edlbAUTO */
90 /**< The indices to communicate, size np */
91 std::vector<gmx_domdec_ind_t> ind;
92 /**< Can we receive data in place? */
93 bool receiveInPlace = false;
96 /*! \brief Load balancing data along a dim used on the master rank of that dim */
101 /**< State var.: max lower bound., incl. neighbors */
102 real cellFracLowerMax = 0;
103 /**< State var.: min upper bound., incl. neighbors */
104 real cellFracUpperMin = 0;
105 /**< Temp. var.: lower limit for cell boundary */
107 /**< Temp. var.: upper limit for cell boundary */
111 /**< Temp. var.: is this cell size at the limit */
112 std::vector<bool> isCellMin;
113 /**< State var.: cell boundaries, box relative */
114 std::vector<real> cellFrac;
115 /**< Temp. var.: old cell size */
116 std::vector<real> oldCellFrac;
118 std::vector<Bounds> bounds;
119 /**< State var.: is DLB limited in this row */
120 bool dlbIsLimited = false;
122 std::vector<real> buf_ncd;
125 /*! \brief Struct for managing cell sizes with DLB along a dimension */
126 struct DDCellsizesWithDlb
128 /**< Cell row root struct, only available on the first rank in a row */
129 std::unique_ptr<RowMaster> rowMaster;
130 /**< The cell sizes, in fractions, along a row, not available on the first rank in a row */
131 std::vector<real> fracRow;
132 /**< The lower corner, in fractions, in triclinic space */
134 /**< The upper corner, in fractions, in triclinic space */
136 /**< The maximum lower corner among all our neighbors */
137 real fracLowerMax = 0;
138 /**< The minimum upper corner among all our neighbors */
139 real fracUpperMin = 0;
142 /*! \brief Struct for compute load commuication
144 * Here floats are accurate enough, since these variables
145 * only influence the load balancing, not the actual MD results.
149 /**< The number of load recordings */
151 /**< Scan of the sum of load over dimensions */
152 float* load = nullptr;
153 /**< The sum of the load over the ranks up to our current dimension */
155 /**< The maximum over the ranks contributing to \p sum */
157 /**< Like \p sum, but takes the maximum when the load balancing is limited */
159 /**< Minimum cell volume, relative to the box */
161 /**< The PP time during which PME can overlap */
163 /**< The PME-only rank load */
165 /**< Bit flags that tell if DLB was limited, per dimension */
169 /*! \brief Data needed to sort an atom to the desired location in the local state */
172 /**< Neighborsearch grid cell index */
174 /**< Global atom/charge group index */
176 /**< Local atom/charge group index */
180 /*! \brief Temporary buffers for sorting atoms */
183 /**< Sorted array of indices */
184 std::vector<gmx_cgsort_t> sorted;
185 /**< Array of stationary atom/charge group indices */
186 std::vector<gmx_cgsort_t> stationary;
187 /**< Array of moved atom/charge group indices */
188 std::vector<gmx_cgsort_t> moved;
189 /**< Integer buffer for sorting */
190 std::vector<int> intBuffer;
193 /*! \brief Manages atom ranges and order for the local state atom vectors */
197 /*! \brief The local state atom order
199 * This enum determines the order of the atoms in the local state.
200 * ddnatHOME and ddnatZONE should be first and second,
201 * the others can be ordered as wanted.
203 enum class Type : int
205 Home, /**< The home atoms */
206 Zones, /**< All zones in the eighth shell */
207 Vsites, /**< Atoms communicated for virtual sites */
208 Constraints, /**< Atoms communicated for constraints */
209 Number /**< Not a count, only present for convenience */
212 /*! \brief Returns the start atom index for range \p rangeType */
213 int start(Type rangeType) const
215 if (rangeType == Type::Home)
221 return end_[static_cast<int>(rangeType) - 1];
225 /*! \brief Returns the end atom index for range \p rangeType */
226 int end(Type rangeType) const { return end_[static_cast<int>(rangeType)]; }
228 /*! \brief Returns the number of home atoms */
229 int numHomeAtoms() const { return end_[static_cast<int>(Type::Home)]; }
231 /*! \brief Returns the total number of atoms */
232 int numAtomsTotal() const { return end_[static_cast<int>(Type::Number) - 1]; }
234 /*! \brief Sets the end index of range \p rangeType to \p end
236 * This should be called either with Type::Home or with a type
237 * that is larger than that passed in the previous call to setEnd.
238 * A release assertion for these conditions are present.
240 void setEnd(Type rangeType, int end)
242 GMX_RELEASE_ASSERT(rangeType == Type::Home || rangeType > lastTypeSet_,
243 "Can only set either home or a larger type than the last one");
245 for (int i = static_cast<int>(rangeType); i < static_cast<int>(Type::Number); i++)
250 lastTypeSet_ = rangeType;
254 /*! \brief The list of end atom indices */
255 std::array<int, static_cast<int>(Type::Number)> end_;
256 Type lastTypeSet_ = Type::Number;
259 /*! \brief Enum of dynamic load balancing states
261 * Allowed DLB states and transitions
262 * - intialization at startup:
263 * -> offUser ("-dlb no")
264 * -> onUser ("-dlb yes")
265 * -> offCanTurnOn ("-dlb auto")
267 * - in automatic mode (i.e. initial state offCanTurnOn):
268 * offCanTurnOn -> onCanTurnOff
269 * offCanTurnOn -> offForever
270 * offCanTurnOn -> offTemporarilyLocked
271 * offTemporarilyLocked -> offCanTurnOn
272 * onCanTurnOff -> offCanTurnOn
276 offUser, /**< DLB is permanently off per user request */
277 offForever, /**< DLB is off due to a runtime condition (not supported or causes performance loss) and will never be turned on */
278 offCanTurnOn, /**< DLB is off and will turn on on imbalance */
279 offTemporarilyLocked, /**< DLB is off and temporarily can't turn on */
280 onCanTurnOff, /**< DLB is on and can turn off when slow */
281 onUser, /**< DLB is permanently on per user request */
282 nr /**< The number of DLB states */
285 /*! \brief The PME domain decomposition for one dimension */
288 /**< The dimension */
290 /**< Tells if DD and PME dims match */
291 gmx_bool dim_match = false;
292 /**< The number of PME ranks/domains in this dimension */
294 /**< Cell sizes for determining the PME comm. with SLB */
295 real* slb_dim_f = nullptr;
296 /**< The minimum pp node location, size nslab */
297 int* pp_min = nullptr;
298 /**< The maximum pp node location, size nslab */
299 int* pp_max = nullptr;
300 /**< The maximum shift for coordinate redistribution in PME */
306 /**< The minimum bottom of this zone */
308 /**< The maximum top of this zone */
310 /**< The minimum top of this zone */
312 /**< The maximum bottom communicaton height for this zone */
314 /**< The maximum top communicaton height for this zone */
316 /**< The bottom value of the first cell in this zone */
318 /**< The top value of the first cell in this zone */
320 /**< Bool disguised as a real, 1 when the above data has been set. 0 otherwise */
324 /*! \brief The number of reals in gmx_ddzone_t */
325 constexpr int c_ddzoneNumReals = 8;
328 class DDBufferAccess;
330 /*! \brief Temporary storage container that minimizes (re)allocation and clearing
332 * This is only the storage, actual access happens through DDBufferAccess.
333 * All methods check if the buffer is (not) in use.
339 /*! \brief Returns a buffer of size \p numElements, the elements are undefined */
340 gmx::ArrayRef<T> resize(size_t numElements)
342 GMX_ASSERT(isInUse_, "Should only operate on acquired buffers");
344 if (numElements > buffer_.size())
346 buffer_.resize(numElements);
349 return gmx::arrayRefFromArray(buffer_.data(), numElements);
352 /*! \brief Acquire the buffer for use with size set to \p numElements, the elements are undefined */
353 gmx::ArrayRef<T> acquire(size_t numElements)
355 GMX_RELEASE_ASSERT(!isInUse_, "Should only request free buffers");
358 return resize(numElements);
361 /*! \brief Releases the buffer, buffer_ should not be used after this */
364 GMX_RELEASE_ASSERT(isInUse_, "Should only release buffers in use");
368 std::vector<T> buffer_; /**< The actual memory buffer */
369 bool isInUse_ = false; /**< Flag that tells whether the buffer is in use */
371 friend class DDBufferAccess<T>;
374 /*! \brief Class that manages access to a temporary memory buffer */
379 /*! \brief Constructor, returns a buffer of size \p numElements, element values are undefined
381 * \note The actual memory buffer \p ddBuffer can not be used to
382 * create other DDBufferAccess objects until the one created
385 DDBufferAccess(DDBuffer<T>& ddBuffer, size_t numElements) : ddBuffer_(ddBuffer)
387 buffer = ddBuffer_.acquire(numElements);
390 ~DDBufferAccess() { ddBuffer_.release(); }
392 /*! \brief Resizes the buffer to \p numElements, new elements are undefined
394 * \note The buffer arrayref is updated after this call.
396 void resize(size_t numElements) { buffer = ddBuffer_.resize(numElements); }
399 DDBuffer<T>& ddBuffer_; /**< Reference to the storage class */
401 gmx::ArrayRef<T> buffer; /**< The access to the memory buffer */
404 /*! \brief Temporary buffer for setting up communiation over one pulse and all zones in the halo */
405 struct dd_comm_setup_work_t
407 /**< The local atom group indices to send */
408 std::vector<int> localAtomGroupBuffer;
409 /**< Buffer for collecting the global atom group indices to send */
410 std::vector<int> atomGroupBuffer;
411 /**< Buffer for collecting the atom group positions to send */
412 std::vector<gmx::RVec> positionBuffer;
413 /**< The number of atoms contained in the atom groups to send */
415 /**< The number of atom groups to send for the last zone */
419 /*! \brief Information about the simulated system */
422 //! True when update groups are used
423 bool useUpdateGroups = false;
424 //! Update atom grouping for each molecule type
425 std::vector<gmx::RangePartitioning> updateGroupingPerMoleculetype;
426 //! The maximum radius over all update groups
427 real maxUpdateGroupRadius;
429 //! Are molecules always whole, i.e. not broken by PBC?
430 bool moleculesAreAlwaysWhole = false;
431 //! Are there inter-domain bonded interactions?
432 bool haveInterDomainBondeds = false;
433 //! Are there inter-domain multi-body interactions?
434 bool haveInterDomainMultiBodyBondeds = false;
436 //! Cut-off for multi-body interactions
437 real minCutoffForMultiBody = 0;
438 //! Cut-off for non-bonded/2-body interactions
440 //! The lower limit for the DD cell size
441 real cellsizeLimit = 0;
443 //! Can atoms connected by constraints be assigned to different domains?
444 bool haveSplitConstraints = false;
445 //! Can atoms connected by settles be assigned to different domains?
446 bool haveSplitSettles = false;
447 //! Estimated communication range needed for constraints
448 real constraintCommunicationRange = 0;
450 //! Whether to only communicate atoms beyond the non-bonded cut-off when they are involved in bonded interactions with non-local atoms
451 bool filterBondedCommunication = false;
452 //! Whether to increase the multi-body cut-off beyond the minimum required
453 bool increaseMultiBodyCutoff = false;
456 /*! \brief Settings that affect the behavior of the domain decomposition
458 * These settings depend on options chosen by the user, set by enviroment
459 * variables, as well as hardware support. The initial DLB state also
460 * depends on the integrator.
462 * Note: Settings that depend on the simulated system are in DDSystemInfo.
466 //! Use MPI_Sendrecv communication instead of non-blocking calls
467 bool useSendRecv2 = false;
469 /* Information for managing the dynamic load balancing */
470 //! Maximum DLB scaling per load balancing step in percent
471 int dlb_scale_lim = 0;
472 //! Flop counter (0=no,1=yes,2=with (eFlop-1)*5% noise
475 //! Request 1D domain decomposition with maximum one communication pulse
476 bool request1DAnd1Pulse;
478 //! Whether to order the DD dimensions from z to x
479 bool useDDOrderZYX = false;
481 //! Whether to use MPI Cartesian reordering of communicators, when supported (almost never)
482 bool useCartesianReorder = true;
484 //! Whether we should record the load
485 bool recordLoad = false;
488 //! Step interval for dumping the local+non-local atoms to pdb
490 //! Step interval for duming the DD grid to pdb
491 int nstDDDumpGrid = 0;
492 //! DD debug print level: 0, 1, 2
495 //! The DLB state at the start of the run
496 DlbState initialDlbState = DlbState::offCanTurnOn;
499 /*! \brief Information on how the DD ranks are set up */
502 /**< The number of particle-particle (non PME-only) ranks */
504 /**< The DD PP grid */
505 ivec numPPCells = { 0, 0, 0 };
507 /* PME and Cartesian communicator stuff */
508 bool usePmeOnlyRanks = false;
509 /**< The number of decomposition dimensions for PME, 0: no PME */
510 int npmedecompdim = 0;
511 /**< The number of ranks doing PME (PP/PME or only PME) */
512 int numRanksDoingPme = 0;
513 /**< The number of PME ranks/domains along x */
515 /**< The number of PME ranks/domains along y */
517 /**< The 1D or 2D PME domain decomposition setup */
518 gmx_ddpme_t ddpme[2];
521 /*! \brief Information on Cartesian MPI setup of the DD ranks */
522 struct CartesianRankSetup
524 /**< Use Cartesian communication between PP and PME ranks */
525 bool bCartesianPP_PME = false;
526 /**< Cartesian grid for combinted PP+PME ranks */
528 /**< The number of dimensions for the PME setup that are Cartesian */
530 /**< The Cartesian index to sim rank conversion, used with bCartesianPP_PME */
531 std::vector<int> ddindex2simnodeid;
533 /* The DD particle-particle nodes only */
534 /**< Use a Cartesian communicator for PP */
535 bool bCartesianPP = false;
536 /**< The Cartesian index to DD rank conversion, used with bCartesianPP */
537 std::vector<int> ddindex2ddnodeid;
540 /*! \brief Struct for domain decomposition communication
542 * This struct contains most information about domain decomposition
543 * communication setup, some communication buffers, some statistics
544 * and also the setup for the communication between particle-particle
545 * and PME only ranks.
547 * All arrays are indexed with 0 to dd->ndim (not Cartesian indexing),
548 * unless stated otherwise.
550 struct gmx_domdec_comm_t // NOLINT (clang-analyzer-optin.performance.Padding)
552 /**< Constant parameters that control DD behavior */
553 DDSettings ddSettings;
555 /**< Information on how the DD ranks are set up */
556 DDRankSetup ddRankSetup;
557 /**< Information on the Cartesian part of the DD rank setup */
558 CartesianRankSetup cartesianRankSetup;
560 /* Charge group / atom sorting */
561 /**< Data structure for cg/atom sorting */
562 std::unique_ptr<gmx_domdec_sort_t> sort;
564 //! Centers of mass of local update groups
565 std::unique_ptr<gmx::UpdateGroupsCog> updateGroupsCog;
567 /* Data for the optional filtering of communication of atoms for bonded interactions */
568 /**< Links between atoms through bonded interactions */
569 t_blocka* bondedLinks = nullptr;
571 /* The DLB state, possible values are defined above */
573 /* With dlbState=DlbState::offCanTurnOn, should we check if to DLB on at the next DD? */
574 gmx_bool bCheckWhetherToTurnDlbOn = false;
575 /* The first DD count since we are running without DLB */
576 int ddPartioningCountFirstDlbOff = 0;
578 /* Cell sizes for static load balancing, first index cartesian */
579 real** slb_frac = nullptr;
581 /**< Information about the simulated system */
582 DDSystemInfo systemInfo;
584 /* The width of the communicated boundaries */
585 /**< Cut-off for multi-body interactions, also 2-body bonded when \p cutoff_mody > \p cutoff */
586 real cutoff_mbody = 0;
587 /**< The minimum guaranteed cell-size, Cartesian indexing */
588 rvec cellsize_min = {};
589 /**< The minimum guaranteed cell-size with dlb=auto */
590 rvec cellsize_min_dlb = {};
591 /**< The lower limit for the DD cell size with DLB */
592 real cellsize_limit = 0;
593 /**< Effectively no NB cut-off limit with DLB for systems without PBC? */
594 gmx_bool bVacDLBNoLimit = false;
596 /** With PME load balancing we set limits on DLB */
597 gmx_bool bPMELoadBalDLBLimits = false;
598 /** DLB needs to take into account that we want to allow this maximum
599 * cut-off (for PME load balancing), this could limit cell boundaries.
601 real PMELoadBal_max_cutoff = 0;
603 /**< box lower corner, required with dim's without pbc and -gcom */
605 /**< box size, required with dim's without pbc and -gcom */
608 /**< The DD cell lower corner, in triclinic space */
610 /**< The DD cell upper corner, in triclinic space */
613 /**< The old \p cell_x0, to check cg displacements */
614 rvec old_cell_x0 = {};
615 /**< The old \p cell_x1, to check cg displacements */
616 rvec old_cell_x1 = {};
618 /** The communication setup and charge group boundaries for the zones */
619 gmx_domdec_zones_t zones;
621 /* The zone limits for DD dimensions 1 and 2 (not 0), determined from
622 * cell boundaries of neighboring cells for staggered grids when using
623 * dynamic load balancing.
625 /**< Zone limits for dim 1 with staggered grids */
626 gmx_ddzone_t zone_d1[2];
627 /**< Zone limits for dim 2 with staggered grids */
628 gmx_ddzone_t zone_d2[2][2];
630 /** The coordinate/force communication setup and indices */
631 gmx_domdec_comm_dim_t cd[DIM];
632 /** Restricts the maximum number of cells to communicate with in one dimension
634 * Dynamic load balancing is not permitted to change sizes if it
635 * would violate this restriction. */
638 /** Which cg distribution is stored on the master node,
639 * stored as DD partitioning call count.
641 int64_t master_cg_ddp_count = 0;
643 /** The number of cg's received from the direct neighbors */
644 int zone_ncg1[DD_MAXZONE] = { 0 };
646 /** The atom ranges in the local state */
647 DDAtomRanges atomRanges;
649 /** Array for signalling if atoms have moved to another domain */
650 std::vector<int> movedBuffer;
652 /** Communication int buffer for general use */
653 DDBuffer<int> intBuffer;
655 /** Communication rvec buffer for general use */
656 DDBuffer<gmx::RVec> rvecBuffer;
658 /* Temporary storage for thread parallel communication setup */
659 /**< Thread-local work data */
660 std::vector<dd_comm_setup_work_t> dth;
662 /* Communication buffer only used with multiple grid pulses */
663 /**< Another rvec comm. buffer */
664 DDBuffer<gmx::RVec> rvecBuffer2;
666 /* Communication buffers for local redistribution */
667 /**< Charge group flag comm. buffers */
668 std::array<std::vector<int>, DIM * 2> cggl_flag;
669 /**< Charge group center comm. buffers */
670 std::array<std::vector<gmx::RVec>, DIM * 2> cgcm_state;
672 /* Cell sizes for dynamic load balancing */
673 std::vector<DDCellsizesWithDlb> cellsizesWithDlb;
675 /* Stuff for load communication */
676 /**< The recorded load data */
677 domdec_load_t* load = nullptr;
678 /**< The number of MPI ranks sharing the GPU our rank is using */
679 int nrank_gpu_shared = 0;
681 /**< The MPI load communicator */
682 MPI_Comm* mpi_comm_load = nullptr;
683 /**< The MPI load communicator for ranks sharing a GPU */
684 MPI_Comm mpi_comm_gpu_shared;
687 /**< Struct for timing the force load balancing region */
688 BalanceRegion* balanceRegion = nullptr;
690 /* Cycle counters over nstlist steps */
691 /**< Total cycles counted */
692 float cycl[ddCyclNr] = {};
693 /**< The number of cycle recordings */
694 int cycl_n[ddCyclNr] = {};
695 /**< The maximum cycle count */
696 float cycl_max[ddCyclNr] = {};
697 /**< Total flops counted */
699 /**< The number of flop recordings */
701 /** How many times did we have load measurements */
703 /** How many times have we collected the load measurements */
704 int n_load_collect = 0;
706 /* Cycle count history for DLB checks */
707 /**< The averaged cycles per step over the last nstlist step before turning on DLB */
708 float cyclesPerStepBeforeDLB = 0;
709 /**< The running average of the cycles per step during DLB */
710 float cyclesPerStepDlbExpAverage = 0;
711 /**< Have we turned off DLB (after turning DLB on)? */
712 bool haveTurnedOffDlb = false;
713 /**< The DD step at which we last measured that DLB off was faster than DLB on, 0 if there was no such step */
714 int64_t dlbSlowerPartitioningCount = 0;
716 /* Statistics for atoms */
717 /**< The atoms per range, summed over the steps */
718 double sum_nat[static_cast<int>(DDAtomRanges::Type::Number)] = {};
720 /* Statistics for calls and times */
721 /**< The number of partioning calls */
723 /**< The number of load recordings */
725 /**< Total MD step time */
726 double load_step = 0.0;
727 /**< Total PP force time */
728 double load_sum = 0.0;
729 /**< Max \p load_sum over the ranks */
730 double load_max = 0.0;
731 /**< Was load balancing limited, per DD dim */
733 /**< Total time on PP done during PME overlap time */
734 double load_mdf = 0.0;
735 /**< Total time on our PME-only rank */
736 double load_pme = 0.0;
738 /** The last partition step */
739 int64_t partition_step = INT_MIN;
742 /*! \brief DD zone permutation
744 * Zone permutation from the Cartesian x-major/z-minor order to an order
745 * that leads to consecutive charge groups for neighbor searching.
746 * TODO: It should be possible to remove this now that the group scheme is removed
748 static const int zone_perm[3][4] = { { 0, 0, 0, 0 }, { 1, 0, 0, 0 }, { 3, 0, 1, 2 } };
750 /*! \brief DD zone reordering to Cartesian order
752 * Index to reorder the zone such that the end up in Cartesian order
753 * with dimension index 0 major and dimension index 2 minor.
755 static const int zone_reorder_cartesian[DD_MAXZONE] = { 0, 1, 3, 2, 5, 4, 6, 7 };
757 /* dd_zo and dd_zp3 is set up such that i zones with non-zero
758 * components see only j zones with that component 0.
761 /*! \brief Returns the DD cut-off distance for multi-body interactions */
762 real dd_cutoff_multibody(const gmx_domdec_t* dd);
764 /*! \brief Returns the DD cut-off distance for two-body interactions */
765 real dd_cutoff_twobody(const gmx_domdec_t* dd);
767 /*! \brief Returns the domain index given the number of domains and the domain coordinates
769 * This order is required to minimize the coordinate communication in PME
770 * which uses decomposition in the x direction.
772 static inline int dd_index(const ivec numDomains, const ivec domainCoordinates)
774 return ((domainCoordinates[XX] * numDomains[YY] + domainCoordinates[YY]) * numDomains[ZZ])
775 + domainCoordinates[ZZ];
778 /*! Returns the size of the buffer to hold fractional cell boundaries for DD dimension index dimIndex */
779 static inline int ddCellFractionBufferSize(const gmx_domdec_t* dd, int dimIndex)
781 return dd->numCells[dd->dim[dimIndex]] + 1 + dimIndex * 2 + 1 + dimIndex;
784 /*! \brief Maximum number of ranks for using send/recv for state scattering and gathering
786 * Use separate MPI send and receive commands
787 * when #ranks <= c_maxNumRanksUseSendRecvForScatterAndGather
788 * This saves memory (and some copying for small #ranks).
789 * For high parallelization scatter and gather calls are used.
791 static constexpr int c_maxNumRanksUseSendRecvForScatterAndGather = 4;
793 /*! \brief Make DD cells larger by this factor than the limit to avoid rounding issues */
794 static constexpr double DD_CELL_MARGIN = 1.0001;
796 /*! \brief Factor for checking DD cell size limitation during DLB, should be in between 1 and DD_CELL_MARGIN */
797 static constexpr double DD_CELL_MARGIN2 = 1.00005;
799 /*! \brief With pressure scaling, keep cell sizes 2% above the limit to allow for some scaling */
800 static constexpr double DD_PRES_SCALE_MARGIN = 1.02;