1 /*
2 * This file is part of the GROMACS molecular simulation package.
3 *
4 * Copyright (c) 2012,2013,2014,2015,2016 by the GROMACS development team.
5 * Copyright (c) 2017,2018,2019,2020, 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 /*! \libinternal \file
37 * \brief
38 * Functionality for per-atom data in the nbnxm module
39 *
40 * \author Berk Hess <hess@kth.se>
41 * \ingroup module_nbnxm
42 * \inlibraryapi
43 */
44
45
46 #ifndef GMX_NBNXN_ATOMDATA_H
47 #define GMX_NBNXN_ATOMDATA_H
48
49 #include <cstdio>
50
51 #include "gromacs/gpu_utils/devicebuffer_datatype.h"
52 #include "gromacs/gpu_utils/hostallocator.h"
53 #include "gromacs/math/vectypes.h"
54 #include "gromacs/mdtypes/locality.h"
55 #include "gromacs/utility/basedefinitions.h"
56 #include "gromacs/utility/bitmask.h"
57 #include "gromacs/utility/real.h"
58
59 namespace gmx
60 {
61 class MDLogger;
62 }
63
64 struct NbnxmGpu;
65 struct nbnxn_atomdata_t;
66 struct nonbonded_verlet_t;
67 struct tMPI_Atomic;
68
69 class GpuEventSynchronizer;
70
71 namespace Nbnxm
72 {
73 class GridSet;
74 enum class KernelType;
75 } // namespace Nbnxm
76
77 //! Convenience type for vector with aligned memory
78 template<typename T>
79 using AlignedVector = std::vector<T, gmx::AlignedAllocator<T>>;
80
81 enum
82 {
83 nbatXYZ,
84 nbatXYZQ,
85 nbatX4,
86 nbatX8
87 };
88
89 //! Stride for coordinate/force arrays with xyz coordinate storage
90 static constexpr int STRIDE_XYZ = 3;
91 //! Stride for coordinate/force arrays with xyzq coordinate storage
92 static constexpr int STRIDE_XYZQ = 4;
93 //! Size of packs of x, y or z with SIMD 4-grouped packed coordinates/forces
94 static constexpr int c_packX4 = 4;
95 //! Size of packs of x, y or z with SIMD 8-grouped packed coordinates/forces
96 static constexpr int c_packX8 = 8;
97 //! Stridefor a pack of 4 coordinates/forces
98 static constexpr int STRIDE_P4 = DIM * c_packX4;
99 //! Stridefor a pack of 8 coordinates/forces
100 static constexpr int STRIDE_P8 = DIM * c_packX8;
101
102 //! Returns the index in a coordinate array corresponding to atom a
103 template<int packSize>
atom_to_x_index(int a)104 static inline int atom_to_x_index(int a)
105 {
106 return DIM * (a & ~(packSize - 1)) + (a & (packSize - 1));
107 }
108
109 /*! \internal
110 * \brief Struct that holds force and energy output buffers */
111 struct nbnxn_atomdata_output_t
112 {
113 /*! \brief Constructor
114 *
115 * \param[in] kernelType Type of non-bonded kernel
116 * \param[in] numEnergyGroups The number of energy groups
117 * \param[in] simdEnergyBufferStride Stride for entries in the energy buffers for SIMD kernels
118 * \param[in] pinningPolicy Sets the pinning policy for all buffers used on the GPU
119 */
120 nbnxn_atomdata_output_t(Nbnxm::KernelType kernelType,
121 int numEnergyGroups,
122 int simdEnergyBufferStride,
123 gmx::PinningPolicy pinningPolicy);
124
125 //! f, size natoms*fstride
126 gmx::HostVector<real> f;
127 //! Shift force array, size SHIFTS*DIM
128 gmx::HostVector<real> fshift;
129 //! Temporary Van der Waals group energy storage
130 gmx::HostVector<real> Vvdw;
131 //! Temporary Coulomb group energy storage
132 gmx::HostVector<real> Vc;
133 //! Temporary SIMD Van der Waals group energy storage
134 AlignedVector<real> VSvdw;
135 //! Temporary SIMD Coulomb group energy storage
136 AlignedVector<real> VSc;
137 };
138
139 /*! \brief Block size in atoms for the non-bonded thread force-buffer reduction.
140 *
141 * Should be a multiple of all cell and x86 SIMD sizes (i.e. 2, 4 and 8).
142 * Should be small to reduce the reduction and zeroing cost,
143 * but too small will result in overhead.
144 * Currently the block size is NBNXN_BUFFERFLAG_SIZE*3*sizeof(real)=192 bytes.
145 */
146 #if GMX_DOUBLE
147 # define NBNXN_BUFFERFLAG_SIZE 8
148 #else
149 # define NBNXN_BUFFERFLAG_SIZE 16
150 #endif
151
152 /*! \brief We store the reduction flags as gmx_bitmask_t.
153 * This limits the number of flags to BITMASK_SIZE.
154 */
155 #define NBNXN_BUFFERFLAG_MAX_THREADS (BITMASK_SIZE)
156
157
158 /*! \brief LJ combination rules: geometric, Lorentz-Berthelot, none */
159 enum
160 {
161 ljcrGEOM,
162 ljcrLB,
163 ljcrNONE,
164 ljcrNR
165 };
166
167 /*! \internal
168 * \brief Struct that stores atom related data for the nbnxn module
169 *
170 * Note: performance would improve slightly when all std::vector containers
171 * in this struct would not initialize during resize().
172 */
173 struct nbnxn_atomdata_t
174 { //NOLINT(clang-analyzer-optin.performance.Padding)
175 /*! \internal
176 * \brief The actual atom data parameter values */
177 struct Params
178 {
179 /*! \brief Constructor
180 *
181 * \param[in] pinningPolicy Sets the pinning policy for all data that might be transfered to a GPU
182 */
183 Params(gmx::PinningPolicy pinningPolicy);
184
185 //! The number of different atom types
186 int numTypes;
187 //! Lennard-Jone 6*C6 and 12*C12 parameters, size numTypes*2*2
188 gmx::HostVector<real> nbfp;
189 //! Combination rule, see enum defined above
190 int comb_rule;
191 //! LJ parameters per atom type, size numTypes*2
192 gmx::HostVector<real> nbfp_comb;
193 //! As nbfp, but with a stride for the present SIMD architecture
194 AlignedVector<real> nbfp_aligned;
195 //! Atom types per atom
196 gmx::HostVector<int> type;
197 //! LJ parameters per atom for fast SIMD loading
198 gmx::HostVector<real> lj_comb;
199 //! Charges per atom, not set with format nbatXYZQ
200 gmx::HostVector<real> q;
201 //! The number of energy groups
202 int nenergrp;
203 //! 2log(nenergrp)
204 int neg_2log;
205 //! The energy groups, one int entry per cluster, only set when needed
206 gmx::HostVector<int> energrp;
207 };
208
209 /*! \internal
210 * \brief Diagonal and topology exclusion helper data for all SIMD kernels. */
211 struct SimdMasks
212 {
213 SimdMasks();
214
215 //! Helper data for setting up diagonal exclusion masks in the SIMD 4xN kernels
216 AlignedVector<real> diagonal_4xn_j_minus_i;
217 //! Helper data for setting up diaginal exclusion masks in the SIMD 2xNN kernels
218 AlignedVector<real> diagonal_2xnn_j_minus_i;
219 //! Filters for topology exclusion masks for the SIMD kernels
220 AlignedVector<uint32_t> exclusion_filter;
221 //! Filters for topology exclusion masks for double SIMD kernels without SIMD int32 logical support
222 AlignedVector<uint64_t> exclusion_filter64;
223 //! Array of masks needed for exclusions
224 AlignedVector<real> interaction_array;
225 };
226
227 /*! \brief Constructor
228 *
229 * \param[in] pinningPolicy Sets the pinning policy for all data that might be transfered to a GPU
230 */
231 nbnxn_atomdata_t(gmx::PinningPolicy pinningPolicy);
232
233 //! Returns a const reference to the parameters
paramsnbnxn_atomdata_t234 const Params& params() const { return params_; }
235
236 //! Returns a non-const reference to the parameters
paramsDeprecatednbnxn_atomdata_t237 Params& paramsDeprecated() { return params_; }
238
239 //! Returns the current total number of atoms stored
numAtomsnbnxn_atomdata_t240 int numAtoms() const { return numAtoms_; }
241
242 //! Return the coordinate buffer, and q with xFormat==nbatXYZQ
xnbnxn_atomdata_t243 gmx::ArrayRef<const real> x() const { return x_; }
244
245 //! Return the coordinate buffer, and q with xFormat==nbatXYZQ
xnbnxn_atomdata_t246 gmx::ArrayRef<real> x() { return x_; }
247
248 //! Resizes the coordinate buffer and sets the number of atoms
249 void resizeCoordinateBuffer(int numAtoms);
250
251 //! Resizes the force buffers for the current number of atoms
252 void resizeForceBuffers();
253
254 private:
255 //! The LJ and charge parameters
256 Params params_;
257 //! The total number of atoms currently stored
258 int numAtoms_;
259
260 public:
261 //! Number of local atoms
262 int natoms_local;
263 //! The format of x (and q), enum
264 int XFormat;
265 //! The format of f, enum
266 int FFormat;
267 //! Do we need to update shift_vec every step?
268 gmx_bool bDynamicBox;
269 //! Shift vectors, copied from t_forcerec
270 gmx::HostVector<gmx::RVec> shift_vec;
271 //! stride for a coordinate in x (usually 3 or 4)
272 int xstride;
273 //! stride for a coordinate in f (usually 3 or 4)
274 int fstride;
275
276 private:
277 //! x and possibly q, size natoms*xstride
278 gmx::HostVector<real> x_;
279
280 public:
281 //! Masks for handling exclusions in the SIMD kernels
282 const SimdMasks simdMasks;
283
284 //! Output data structures, 1 per thread
285 std::vector<nbnxn_atomdata_output_t> out;
286
287 //! Reduction related data
288 //! \{
289 //! Use the flags or operate on all atoms
290 gmx_bool bUseBufferFlags;
291 //! Flags for buffer zeroing+reduc.
292 std::vector<gmx_bitmask_t> buffer_flags;
293 //! Use tree for force reduction
294 gmx_bool bUseTreeReduce;
295 //! Synchronization step for tree reduce
296 tMPI_Atomic* syncStep;
297 //! \}
298 };
299
300 /*! \brief Copy na rvec elements from x to xnb using nbatFormat, start dest a0,
301 * and fills up to na_round with coordinates that are far away.
302 */
303 void copy_rvec_to_nbat_real(const int* a, int na, int na_round, const rvec* x, int nbatFormat, real* xnb, int a0);
304
305 //! Describes the combination rule in use by this force field
306 enum
307 {
308 enbnxninitcombruleDETECT,
309 enbnxninitcombruleGEOM,
310 enbnxninitcombruleLB,
311 enbnxninitcombruleNONE
312 };
313
314 /*! \brief Initialize the non-bonded atom data structure.
315 *
316 * The enum for nbatXFormat is in the file defining nbnxn_atomdata_t.
317 * Copy the ntypes*ntypes*2 sized nbfp non-bonded parameter list
318 * to the atom data structure.
319 * enbnxninitcombrule sets what combination rule data gets stored in nbat.
320 */
321 void nbnxn_atomdata_init(const gmx::MDLogger& mdlog,
322 nbnxn_atomdata_t* nbat,
323 Nbnxm::KernelType kernelType,
324 int enbnxninitcombrule,
325 int ntype,
326 gmx::ArrayRef<const real> nbfp,
327 int n_energygroups,
328 int nout);
329
330 //! Sets the atomdata after pair search
331 void nbnxn_atomdata_set(nbnxn_atomdata_t* nbat,
332 const Nbnxm::GridSet& gridSet,
333 gmx::ArrayRef<const int> atomTypes,
334 gmx::ArrayRef<const real> atomCharges,
335 gmx::ArrayRef<const int> atomInfo);
336
337 //! Copy the shift vectors to nbat
338 void nbnxn_atomdata_copy_shiftvec(gmx_bool dynamic_box, rvec* shift_vec, nbnxn_atomdata_t* nbat);
339
340 /*! \brief Transform coordinates to xbat layout
341 *
342 * Creates a copy of the coordinates buffer using short-range ordering.
343 *
344 * \param[in] gridSet The grids data.
345 * \param[in] locality If the transformation should be applied to local or non local coordinates.
346 * \param[in] fillLocal Tells if the local filler particle coordinates should be zeroed.
347 * \param[in] coordinates Coordinates in plain rvec format.
348 * \param[in,out] nbat Data in NBNXM format, used for mapping formats and to locate the output buffer.
349 */
350 void nbnxn_atomdata_copy_x_to_nbat_x(const Nbnxm::GridSet& gridSet,
351 gmx::AtomLocality locality,
352 bool fillLocal,
353 const rvec* coordinates,
354 nbnxn_atomdata_t* nbat);
355
356 /*! \brief Transform coordinates to xbat layout on GPU
357 *
358 * Creates a GPU copy of the coordinates buffer using short-range ordering.
359 * As input, uses coordinates in plain rvec format in GPU memory.
360 *
361 * \param[in] gridSet The grids data.
362 * \param[in] locality If the transformation should be applied to local or non local coordinates.
363 * \param[in] fillLocal Tells if the local filler particle coordinates should be zeroed.
364 * \param[in,out] gpu_nbv The NBNXM GPU data structure.
365 * \param[in] d_x Coordinates to be copied (in plain rvec format).
366 * \param[in] xReadyOnDevice Event synchronizer indicating that the coordinates are ready in the device memory.
367 */
368 void nbnxn_atomdata_x_to_nbat_x_gpu(const Nbnxm::GridSet& gridSet,
369 gmx::AtomLocality locality,
370 bool fillLocal,
371 NbnxmGpu* gpu_nbv,
372 DeviceBuffer<gmx::RVec> d_x,
373 GpuEventSynchronizer* xReadyOnDevice);
374
375 /*! \brief Add the computed forces to \p f, an internal reduction might be performed as well
376 *
377 * \param[in] nbat Atom data in NBNXM format.
378 * \param[in] locality If the reduction should be performed on local or non-local atoms.
379 * \param[in] gridSet The grids data.
380 * \param[out] totalForce Buffer to accumulate resulting force
381 */
382 void reduceForces(nbnxn_atomdata_t* nbat, gmx::AtomLocality locality, const Nbnxm::GridSet& gridSet, rvec* totalForce);
383
384 //! Add the fshift force stored in nbat to fshift
385 void nbnxn_atomdata_add_nbat_fshift_to_fshift(const nbnxn_atomdata_t& nbat, gmx::ArrayRef<gmx::RVec> fshift);
386
387 //! Get the atom start index and number of atoms for a given locality
388 void nbnxn_get_atom_range(gmx::AtomLocality atomLocality,
389 const Nbnxm::GridSet& gridSet,
390 int* atomStart,
391 int* nAtoms);
392 #endif
393