2 * This file is part of the GROMACS molecular simulation package.
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.
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.
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.
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.
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.
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.
39 * Declares inline-friendly code for making 2xNN pairlists
41 * \author Berk Hess <hess@kth.se>
42 * \ingroup module_nbnxm
46 //! Stride of the packed x coordinate array
47 static constexpr int c_xStride2xNN
= (GMX_SIMD_REAL_WIDTH
>= 2 * c_nbnxnCpuIClusterSize
)
48 ? GMX_SIMD_REAL_WIDTH
/ 2
49 : c_nbnxnCpuIClusterSize
;
51 //! Copies PBC shifted i-cell packed atom coordinates to working array
52 static inline void icell_set_x_simd_2xnn(int ci
,
56 int gmx_unused stride
,
58 NbnxnPairlistCpuWork
* work
)
61 real
* x_ci_simd
= work
->iClusterData
.xSimd
.data();
63 ia
= xIndexFromCi
<NbnxnLayout::Simd2xNN
>(ci
);
65 store(x_ci_simd
+ 0 * GMX_SIMD_REAL_WIDTH
,
66 loadU1DualHsimd(x
+ ia
+ 0 * c_xStride2xNN
+ 0) + SimdReal(shx
));
67 store(x_ci_simd
+ 1 * GMX_SIMD_REAL_WIDTH
,
68 loadU1DualHsimd(x
+ ia
+ 1 * c_xStride2xNN
+ 0) + SimdReal(shy
));
69 store(x_ci_simd
+ 2 * GMX_SIMD_REAL_WIDTH
,
70 loadU1DualHsimd(x
+ ia
+ 2 * c_xStride2xNN
+ 0) + SimdReal(shz
));
71 store(x_ci_simd
+ 3 * GMX_SIMD_REAL_WIDTH
,
72 loadU1DualHsimd(x
+ ia
+ 0 * c_xStride2xNN
+ 2) + SimdReal(shx
));
73 store(x_ci_simd
+ 4 * GMX_SIMD_REAL_WIDTH
,
74 loadU1DualHsimd(x
+ ia
+ 1 * c_xStride2xNN
+ 2) + SimdReal(shy
));
75 store(x_ci_simd
+ 5 * GMX_SIMD_REAL_WIDTH
,
76 loadU1DualHsimd(x
+ ia
+ 2 * c_xStride2xNN
+ 2) + SimdReal(shz
));
79 /*! \brief SIMD code for checking and adding cluster-pairs to the list using coordinates in packed format.
81 * Checks bouding box distances and possibly atom pair distances.
82 * This is an accelerated version of make_cluster_list_simple.
84 * \param[in] jGrid The j-grid
85 * \param[in,out] nbl The pair-list to store the cluster pairs in
86 * \param[in] icluster The index of the i-cluster
87 * \param[in] firstCell The first cluster in the j-range, using i-cluster size indexing
88 * \param[in] lastCell The last cluster in the j-range, using i-cluster size indexing
89 * \param[in] excludeSubDiagonal Exclude atom pairs with i-index > j-index
90 * \param[in] x_j Coordinates for the j-atom, in SIMD packed format
91 * \param[in] rlist2 The squared list cut-off
92 * \param[in] rbb2 The squared cut-off for putting cluster-pairs in the list based on bounding box distance only
93 * \param[in,out] numDistanceChecks The number of distance checks performed
95 static inline void makeClusterListSimd2xnn(const Grid
& jGrid
,
96 NbnxnPairlistCpu
* nbl
,
100 bool excludeSubDiagonal
,
101 const real
* gmx_restrict x_j
,
104 int* gmx_restrict numDistanceChecks
)
107 const real
* gmx_restrict x_ci_simd
= nbl
->work
->iClusterData
.xSimd
.data();
108 const BoundingBox
* gmx_restrict bb_ci
= nbl
->work
->iClusterData
.bb
.data();
110 SimdReal jx_S
, jy_S
, jz_S
;
112 SimdReal dx_S0
, dy_S0
, dz_S0
;
113 SimdReal dx_S2
, dy_S2
, dz_S2
;
128 int jclusterFirst
= cjFromCi
<NbnxnLayout::Simd2xNN
, 0>(firstCell
);
129 int jclusterLast
= cjFromCi
<NbnxnLayout::Simd2xNN
, 1>(lastCell
);
130 GMX_ASSERT(jclusterLast
>= jclusterFirst
,
131 "We should have a non-empty j-cluster range, since the calling code should have "
132 "ensured a non-empty cell range");
134 rc2_S
= SimdReal(rlist2
);
137 while (!InRange
&& jclusterFirst
<= jclusterLast
)
139 d2
= clusterBoundingBoxDistance2(bb_ci
[0], jGrid
.jBoundingBoxes()[jclusterFirst
]);
140 *numDistanceChecks
+= 2;
142 /* Check if the distance is within the distance where
143 * we use only the bounding box distance rbb,
144 * or within the cut-off and there is at least one atom pair
145 * within the cut-off.
151 else if (d2
< rlist2
)
153 xind_f
= xIndexFromCj
<NbnxnLayout::Simd2xNN
>(
154 cjFromCi
<NbnxnLayout::Simd2xNN
, 0>(jGrid
.cellOffset()) + jclusterFirst
);
156 jx_S
= loadDuplicateHsimd(x_j
+ xind_f
+ 0 * c_xStride2xNN
);
157 jy_S
= loadDuplicateHsimd(x_j
+ xind_f
+ 1 * c_xStride2xNN
);
158 jz_S
= loadDuplicateHsimd(x_j
+ xind_f
+ 2 * c_xStride2xNN
);
160 /* Calculate distance */
161 dx_S0
= load
<SimdReal
>(x_ci_simd
+ 0 * GMX_SIMD_REAL_WIDTH
) - jx_S
;
162 dy_S0
= load
<SimdReal
>(x_ci_simd
+ 1 * GMX_SIMD_REAL_WIDTH
) - jy_S
;
163 dz_S0
= load
<SimdReal
>(x_ci_simd
+ 2 * GMX_SIMD_REAL_WIDTH
) - jz_S
;
164 dx_S2
= load
<SimdReal
>(x_ci_simd
+ 3 * GMX_SIMD_REAL_WIDTH
) - jx_S
;
165 dy_S2
= load
<SimdReal
>(x_ci_simd
+ 4 * GMX_SIMD_REAL_WIDTH
) - jy_S
;
166 dz_S2
= load
<SimdReal
>(x_ci_simd
+ 5 * GMX_SIMD_REAL_WIDTH
) - jz_S
;
168 /* rsq = dx*dx+dy*dy+dz*dz */
169 rsq_S0
= norm2(dx_S0
, dy_S0
, dz_S0
);
170 rsq_S2
= norm2(dx_S2
, dy_S2
, dz_S2
);
172 wco_S0
= (rsq_S0
< rc2_S
);
173 wco_S2
= (rsq_S2
< rc2_S
);
175 wco_any_S
= wco_S0
|| wco_S2
;
177 InRange
= anyTrue(wco_any_S
);
179 *numDistanceChecks
+= 2 * GMX_SIMD_REAL_WIDTH
;
192 while (!InRange
&& jclusterLast
> jclusterFirst
)
194 d2
= clusterBoundingBoxDistance2(bb_ci
[0], jGrid
.jBoundingBoxes()[jclusterLast
]);
195 *numDistanceChecks
+= 2;
197 /* Check if the distance is within the distance where
198 * we use only the bounding box distance rbb,
199 * or within the cut-off and there is at least one atom pair
200 * within the cut-off.
206 else if (d2
< rlist2
)
208 xind_l
= xIndexFromCj
<NbnxnLayout::Simd2xNN
>(
209 cjFromCi
<NbnxnLayout::Simd2xNN
, 0>(jGrid
.cellOffset()) + jclusterLast
);
211 jx_S
= loadDuplicateHsimd(x_j
+ xind_l
+ 0 * c_xStride2xNN
);
212 jy_S
= loadDuplicateHsimd(x_j
+ xind_l
+ 1 * c_xStride2xNN
);
213 jz_S
= loadDuplicateHsimd(x_j
+ xind_l
+ 2 * c_xStride2xNN
);
215 /* Calculate distance */
216 dx_S0
= load
<SimdReal
>(x_ci_simd
+ 0 * GMX_SIMD_REAL_WIDTH
) - jx_S
;
217 dy_S0
= load
<SimdReal
>(x_ci_simd
+ 1 * GMX_SIMD_REAL_WIDTH
) - jy_S
;
218 dz_S0
= load
<SimdReal
>(x_ci_simd
+ 2 * GMX_SIMD_REAL_WIDTH
) - jz_S
;
219 dx_S2
= load
<SimdReal
>(x_ci_simd
+ 3 * GMX_SIMD_REAL_WIDTH
) - jx_S
;
220 dy_S2
= load
<SimdReal
>(x_ci_simd
+ 4 * GMX_SIMD_REAL_WIDTH
) - jy_S
;
221 dz_S2
= load
<SimdReal
>(x_ci_simd
+ 5 * GMX_SIMD_REAL_WIDTH
) - jz_S
;
223 /* rsq = dx*dx+dy*dy+dz*dz */
224 rsq_S0
= norm2(dx_S0
, dy_S0
, dz_S0
);
225 rsq_S2
= norm2(dx_S2
, dy_S2
, dz_S2
);
227 wco_S0
= (rsq_S0
< rc2_S
);
228 wco_S2
= (rsq_S2
< rc2_S
);
230 wco_any_S
= wco_S0
|| wco_S2
;
232 InRange
= anyTrue(wco_any_S
);
234 *numDistanceChecks
+= 2 * GMX_SIMD_REAL_WIDTH
;
242 if (jclusterFirst
<= jclusterLast
)
244 for (int jcluster
= jclusterFirst
; jcluster
<= jclusterLast
; jcluster
++)
246 /* Store cj and the interaction mask */
248 cjEntry
.cj
= cjFromCi
<NbnxnLayout::Simd2xNN
, 0>(jGrid
.cellOffset()) + jcluster
;
249 cjEntry
.excl
= get_imask_simd_2xnn(excludeSubDiagonal
, icluster
, jcluster
);
250 nbl
->cj
.push_back(cjEntry
);
252 /* Increase the closing index in the i list */
253 nbl
->ci
.back().cj_ind_end
= nbl
->cj
.size();