2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,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.
36 #ifndef GMX_NBNXM_PAIRLIST_H
37 #define GMX_NBNXM_PAIRLIST_H
43 #include "gromacs/gpu_utils/hostallocator.h"
44 #include "gromacs/math/vectypes.h"
45 #include "gromacs/mdtypes/nblist.h"
46 #include "gromacs/utility/basedefinitions.h"
47 #include "gromacs/utility/defaultinitializationallocator.h"
48 #include "gromacs/utility/enumerationhelpers.h"
49 #include "gromacs/utility/real.h"
51 // This file with constants is separate from this file to be able
52 // to include it during OpenCL jitting without including config.h
53 #include "gromacs/nbnxm/constants.h"
55 struct NbnxnListParameters;
56 struct NbnxnPairlistCpuWork;
57 struct NbnxnPairlistGpuWork;
61 enum class KernelType;
64 /* Convenience type for vector with aligned memory */
66 using AlignedVector = std::vector < T, gmx::AlignedAllocator < T>>;
68 /* Convenience type for vector that avoids initialization at resize() */
70 using FastVector = std::vector < T, gmx::DefaultInitializationAllocator < T>>;
72 enum class PairlistType : int
81 static constexpr gmx::EnumerationArray<PairlistType, int> IClusterSizePerListType = { { 4, 4, 4, 8 } };
82 static constexpr gmx::EnumerationArray<PairlistType, int> JClusterSizePerListType = { { 2, 4, 8, 8 } };
84 /* With CPU kernels the i-cluster size is always 4 atoms. */
85 static constexpr int c_nbnxnCpuIClusterSize = 4;
87 /* With GPU kernels the i and j cluster size is 8 atoms for CUDA and can be set at compile time for OpenCL */
88 #if GMX_GPU == GMX_GPU_OPENCL
89 static constexpr int c_nbnxnGpuClusterSize = GMX_OPENCL_NB_CLUSTER_SIZE;
91 static constexpr int c_nbnxnGpuClusterSize = 8;
94 /* The number of clusters in a pair-search cell, used for GPU */
95 static constexpr int c_gpuNumClusterPerCellZ = 2;
96 static constexpr int c_gpuNumClusterPerCellY = 2;
97 static constexpr int c_gpuNumClusterPerCellX = 2;
98 static constexpr int c_gpuNumClusterPerCell = c_gpuNumClusterPerCellZ*c_gpuNumClusterPerCellY*c_gpuNumClusterPerCellX;
101 /* In CUDA the number of threads in a warp is 32 and we have cluster pairs
102 * of 8*8=64 atoms, so it's convenient to store data for cluster pair halves.
104 static constexpr int c_nbnxnGpuClusterpairSplit = 2;
106 /* The fixed size of the exclusion mask array for a half cluster pair */
107 static constexpr int c_nbnxnGpuExclSize = c_nbnxnGpuClusterSize*c_nbnxnGpuClusterSize/c_nbnxnGpuClusterpairSplit;
109 /* A buffer data structure of 64 bytes
110 * to be placed at the beginning and end of structs
111 * to avoid cache invalidation of the real contents
112 * of the struct by writes to neighboring memory.
116 } gmx_cache_protect_t;
118 /* This is the actual cluster-pair list j-entry.
119 * cj is the j-cluster.
120 * The interaction bits in excl are indexed i-major, j-minor.
121 * The cj entries are sorted such that ones with exclusions come first.
122 * This means that once a full mask (=NBNXN_INTERACTION_MASK_ALL)
123 * is found, all subsequent j-entries in the i-entry also have full masks.
127 int cj; /* The j-cluster */
128 unsigned int excl; /* The exclusion (interaction) bits */
131 /* In nbnxn_ci_t the integer shift contains the shift in the lower 7 bits.
132 * The upper bits contain information for non-bonded kernel optimization.
133 * Simply calculating LJ and Coulomb for all pairs in a cluster pair is fine.
134 * But three flags can be used to skip interactions, currently only for subc=0
135 * !(shift & NBNXN_CI_DO_LJ(subc)) => we can skip LJ for all pairs
136 * shift & NBNXN_CI_HALF_LJ(subc) => we can skip LJ for the second half of i
137 * !(shift & NBNXN_CI_DO_COUL(subc)) => we can skip Coulomb for all pairs
139 #define NBNXN_CI_SHIFT 127
140 #define NBNXN_CI_DO_LJ(subc) (1<<(7+3*(subc)))
141 #define NBNXN_CI_HALF_LJ(subc) (1<<(8+3*(subc)))
142 #define NBNXN_CI_DO_COUL(subc) (1<<(9+3*(subc)))
144 /* Cluster-pair Interaction masks
145 * Bit i*j-cluster-size + j tells if atom i and j interact.
147 // TODO: Rename according to convention when moving into Nbnxn namespace
148 /* All interaction mask is the same for all kernels */
149 constexpr unsigned int NBNXN_INTERACTION_MASK_ALL = 0xffffffffU;
150 /* 4x4 kernel diagonal mask */
151 constexpr unsigned int NBNXN_INTERACTION_MASK_DIAG = 0x08ceU;
152 /* 4x2 kernel diagonal masks */
153 constexpr unsigned int NBNXN_INTERACTION_MASK_DIAG_J2_0 = 0x0002U;
154 constexpr unsigned int NBNXN_INTERACTION_MASK_DIAG_J2_1 = 0x002fU;
155 /* 4x8 kernel diagonal masks */
156 constexpr unsigned int NBNXN_INTERACTION_MASK_DIAG_J8_0 = 0xf0f8fcfeU;
157 constexpr unsigned int NBNXN_INTERACTION_MASK_DIAG_J8_1 = 0x0080c0e0U;
159 /* Simple pair-list i-unit */
162 int ci; /* i-cluster */
163 int shift; /* Shift vector index plus possible flags, see above */
164 int cj_ind_start; /* Start index into cj */
165 int cj_ind_end; /* End index into cj */
168 /* Grouped pair-list i-unit */
170 /* Returns the number of j-cluster groups in this entry */
171 int numJClusterGroups() const
173 return cj4_ind_end - cj4_ind_start;
176 int sci; /* i-super-cluster */
177 int shift; /* Shift vector index plus possible flags */
178 int cj4_ind_start; /* Start index into cj4 */
179 int cj4_ind_end; /* End index into cj4 */
182 /* Interaction data for a j-group for one warp */
185 // The i-cluster interactions mask for 1 warp
186 unsigned int imask = 0U;
187 // Index into the exclusion array for 1 warp, default index 0 which means no exclusions
192 int cj[c_nbnxnGpuJgroupSize]; /* The 4 j-clusters */
193 nbnxn_im_ei_t imei[c_nbnxnGpuClusterpairSplit]; /* The i-cluster mask data for 2 warps */
196 /* Struct for storing the atom-pair interaction bits for a cluster pair in a GPU pairlist */
199 /* Constructor, sets no exclusions, so all atom pairs interacting */
202 for (unsigned int &pairEntry : pair)
204 pairEntry = NBNXN_INTERACTION_MASK_ALL;
208 /* Topology exclusion interaction bits per warp */
209 unsigned int pair[c_nbnxnGpuExclSize];
212 /* Cluster pairlist type for use on CPUs */
213 struct NbnxnPairlistCpu
215 gmx_cache_protect_t cp0;
217 int na_ci; /* The number of atoms per i-cluster */
218 int na_cj; /* The number of atoms per j-cluster */
219 real rlist; /* The radius for constructing the list */
220 FastVector<nbnxn_ci_t> ci; /* The i-cluster list */
221 FastVector<nbnxn_ci_t> ciOuter; /* The outer, unpruned i-cluster list */
223 FastVector<nbnxn_cj_t> cj; /* The j-cluster list, size ncj */
224 FastVector<nbnxn_cj_t> cjOuter; /* The outer, unpruned j-cluster list */
225 int ncjInUse; /* The number of j-clusters that are used by ci entries in this list, will be <= cj.size() */
227 int nci_tot; /* The total number of i clusters */
229 NbnxnPairlistCpuWork *work;
231 gmx_cache_protect_t cp1;
234 /* Cluster pairlist type, with extra hierarchies, for on the GPU
236 * NOTE: for better performance when combining lists over threads,
237 * all vectors should use default initialization. But when
238 * changing this, excl should be intialized when adding entries.
240 struct NbnxnPairlistGpu
244 * \param[in] pinningPolicy Sets the pinning policy for all buffers used on the GPU
246 NbnxnPairlistGpu(gmx::PinningPolicy pinningPolicy);
248 gmx_cache_protect_t cp0;
250 int na_ci; /* The number of atoms per i-cluster */
251 int na_cj; /* The number of atoms per j-cluster */
252 int na_sc; /* The number of atoms per super cluster */
253 real rlist; /* The radius for constructing the list */
254 // The i-super-cluster list, indexes into cj4;
255 gmx::HostVector<nbnxn_sci_t> sci;
256 // The list of 4*j-cluster groups
257 gmx::HostVector<nbnxn_cj4_t> cj4;
258 // Atom interaction bits (non-exclusions)
259 gmx::HostVector<nbnxn_excl_t> excl;
260 // The total number of i-clusters
263 NbnxnPairlistGpuWork *work;
265 gmx_cache_protect_t cp1;
268 struct nbnxn_pairlist_set_t
270 nbnxn_pairlist_set_t(const NbnxnListParameters &listParams);
272 int nnbl; /* number of lists */
273 NbnxnPairlistCpu **nbl; /* lists for CPU */
274 NbnxnPairlistCpu **nbl_work; /* work space for rebalancing lists */
275 NbnxnPairlistGpu **nblGpu; /* lists for GPU */
276 const NbnxnListParameters ¶ms; /* Pairlist parameters desribing setup and ranges */
277 gmx_bool bCombined; /* TRUE if lists get combined into one (the 1st) */
278 gmx_bool bSimple; /* TRUE if the list of of type "simple"
279 (na_sc=na_s, no super-clusters used) */
281 /* Counts for debug printing */
282 int natpair_ljq; /* Total number of atom pairs for LJ+Q kernel */
283 int natpair_lj; /* Total number of atom pairs for LJ kernel */
284 int natpair_q; /* Total number of atom pairs for Q kernel */
285 std::vector<t_nblist *> nbl_fep; /* List of free-energy atom pair interactions */
288 //! Initializes a free-energy pair-list
289 void nbnxn_init_pairlist_fep(t_nblist *nl);