2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014,2015,2016,2017,2018, 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 #if GMX_SIMD_REAL_WIDTH >= NBNXN_CPU_CLUSTER_I_SIZE
37 #define STRIDE_S (GMX_SIMD_REAL_WIDTH)
39 #define STRIDE_S NBNXN_CPU_CLUSTER_I_SIZE
42 /* Copies PBC shifted i-cell packed atom coordinates to working array */
44 icell_set_x_simd_4xn(int ci,
45 real shx, real shy, real shz,
46 int gmx_unused stride, const real *x,
47 nbnxn_list_work_t *work)
50 real *x_ci_simd = work->x_ci_simd;
52 ia = xIndexFromCi<NbnxnLayout::Simd4xN>(ci);
54 store(x_ci_simd + 0*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 0*STRIDE_S ] + shx) );
55 store(x_ci_simd + 1*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 1*STRIDE_S ] + shy) );
56 store(x_ci_simd + 2*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 2*STRIDE_S ] + shz) );
57 store(x_ci_simd + 3*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 0*STRIDE_S + 1] + shx) );
58 store(x_ci_simd + 4*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 1*STRIDE_S + 1] + shy) );
59 store(x_ci_simd + 5*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 2*STRIDE_S + 1] + shz) );
60 store(x_ci_simd + 6*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 0*STRIDE_S + 2] + shx) );
61 store(x_ci_simd + 7*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 1*STRIDE_S + 2] + shy) );
62 store(x_ci_simd + 8*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 2*STRIDE_S + 2] + shz) );
63 store(x_ci_simd + 9*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 0*STRIDE_S + 3] + shx) );
64 store(x_ci_simd + 10*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 1*STRIDE_S + 3] + shy) );
65 store(x_ci_simd + 11*GMX_SIMD_REAL_WIDTH, SimdReal(x[ia + 2*STRIDE_S + 3] + shz) );
68 /* SIMD code for checking and adding cluster-pairs to the list using coordinates in packed format.
70 * Checks bouding box distances and possibly atom pair distances.
71 * This is an accelerated version of make_cluster_list_simple.
73 * \param[in] gridj The j-grid
74 * \param[in,out] nbl The pair-list to store the cluster pairs in
75 * \param[in] icluster The index of the i-cluster
76 * \param[in] firstCell The first cluster in the j-range, using i-cluster size indexing
77 * \param[in] lastCell The last cluster in the j-range, using i-cluster size indexing
78 * \param[in] excludeSubDiagonal Exclude atom pairs with i-index > j-index
79 * \param[in] x_j Coordinates for the j-atom, in SIMD packed format
80 * \param[in] rlist2 The squared list cut-off
81 * \param[in] rbb2 The squared cut-off for putting cluster-pairs in the list based on bounding box distance only
82 * \param[in,out] numDistanceChecks The number of distance checks performed
85 makeClusterListSimd4xn(const nbnxn_grid_t * gridj,
86 nbnxn_pairlist_t * nbl,
90 bool excludeSubDiagonal,
91 const real * gmx_restrict x_j,
94 int * gmx_restrict numDistanceChecks)
97 const real * gmx_restrict x_ci_simd = nbl->work->x_ci_simd;
98 const nbnxn_bb_t * gmx_restrict bb_ci = nbl->work->bb_ci;
100 SimdReal jx_S, jy_S, jz_S;
102 SimdReal dx_S0, dy_S0, dz_S0;
103 SimdReal dx_S1, dy_S1, dz_S1;
104 SimdReal dx_S2, dy_S2, dz_S2;
105 SimdReal dx_S3, dy_S3, dz_S3;
116 SimdBool wco_any_S01, wco_any_S23, wco_any_S;
124 /* Convert the j-range from i-cluster size indexing to j-cluster indexing */
125 /* cppcheck-suppress selfAssignment . selfAssignment for width 4.*/
126 int jclusterFirst = cjFromCi<NbnxnLayout::Simd4xN>(firstCell);
127 #if GMX_SIMD_REAL_WIDTH >= NBNXN_CPU_CLUSTER_I_SIZE
128 int jclusterLast = cjFromCi<NbnxnLayout::Simd4xN>(lastCell);
130 /* Set the correct last j-cluster with a j-cluster size of 2 */
131 int jclusterLast = cjFromCi<NbnxnLayout::Simd4xN>(lastCell + 1) - 1;
133 GMX_ASSERT(jclusterLast >= jclusterFirst, "We should have a non-empty j-cluster range, since the calling code should have ensured a non-empty cell range");
135 rc2_S = SimdReal(rlist2);
138 while (!InRange && jclusterFirst <= jclusterLast)
140 #if NBNXN_SEARCH_BB_SIMD4
141 d2 = subc_bb_dist2_simd4(0, bb_ci, jclusterFirst, gridj->bbj);
143 d2 = subc_bb_dist2(0, bb_ci, jclusterFirst, gridj->bbj);
145 *numDistanceChecks += 2;
147 /* Check if the distance is within the distance where
148 * we use only the bounding box distance rbb,
149 * or within the cut-off and there is at least one atom pair
150 * within the cut-off.
156 else if (d2 < rlist2)
158 xind_f = xIndexFromCj<NbnxnLayout::Simd4xN>(cjFromCi<NbnxnLayout::Simd4xN>(gridj->cell0) + jclusterFirst);
160 jx_S = load<SimdReal>(x_j + xind_f + 0*STRIDE_S);
161 jy_S = load<SimdReal>(x_j + xind_f + 1*STRIDE_S);
162 jz_S = load<SimdReal>(x_j + xind_f + 2*STRIDE_S);
165 /* Calculate distance */
166 dx_S0 = load<SimdReal>(x_ci_simd + 0*GMX_SIMD_REAL_WIDTH) - jx_S;
167 dy_S0 = load<SimdReal>(x_ci_simd + 1*GMX_SIMD_REAL_WIDTH) - jy_S;
168 dz_S0 = load<SimdReal>(x_ci_simd + 2*GMX_SIMD_REAL_WIDTH) - jz_S;
169 dx_S1 = load<SimdReal>(x_ci_simd + 3*GMX_SIMD_REAL_WIDTH) - jx_S;
170 dy_S1 = load<SimdReal>(x_ci_simd + 4*GMX_SIMD_REAL_WIDTH) - jy_S;
171 dz_S1 = load<SimdReal>(x_ci_simd + 5*GMX_SIMD_REAL_WIDTH) - jz_S;
172 dx_S2 = load<SimdReal>(x_ci_simd + 6*GMX_SIMD_REAL_WIDTH) - jx_S;
173 dy_S2 = load<SimdReal>(x_ci_simd + 7*GMX_SIMD_REAL_WIDTH) - jy_S;
174 dz_S2 = load<SimdReal>(x_ci_simd + 8*GMX_SIMD_REAL_WIDTH) - jz_S;
175 dx_S3 = load<SimdReal>(x_ci_simd + 9*GMX_SIMD_REAL_WIDTH) - jx_S;
176 dy_S3 = load<SimdReal>(x_ci_simd + 10*GMX_SIMD_REAL_WIDTH) - jy_S;
177 dz_S3 = load<SimdReal>(x_ci_simd + 11*GMX_SIMD_REAL_WIDTH) - jz_S;
179 /* rsq = dx*dx+dy*dy+dz*dz */
180 rsq_S0 = norm2(dx_S0, dy_S0, dz_S0);
181 rsq_S1 = norm2(dx_S1, dy_S1, dz_S1);
182 rsq_S2 = norm2(dx_S2, dy_S2, dz_S2);
183 rsq_S3 = norm2(dx_S3, dy_S3, dz_S3);
185 wco_S0 = (rsq_S0 < rc2_S);
186 wco_S1 = (rsq_S1 < rc2_S);
187 wco_S2 = (rsq_S2 < rc2_S);
188 wco_S3 = (rsq_S3 < rc2_S);
190 wco_any_S01 = wco_S0 || wco_S1;
191 wco_any_S23 = wco_S2 || wco_S3;
192 wco_any_S = wco_any_S01 || wco_any_S23;
194 InRange = anyTrue(wco_any_S);
196 *numDistanceChecks += 4*GMX_SIMD_REAL_WIDTH;
209 while (!InRange && jclusterLast > jclusterFirst)
211 #if NBNXN_SEARCH_BB_SIMD4
212 d2 = subc_bb_dist2_simd4(0, bb_ci, jclusterLast, gridj->bbj);
214 d2 = subc_bb_dist2(0, bb_ci, jclusterLast, gridj->bbj);
216 *numDistanceChecks += 2;
218 /* Check if the distance is within the distance where
219 * we use only the bounding box distance rbb,
220 * or within the cut-off and there is at least one atom pair
221 * within the cut-off.
227 else if (d2 < rlist2)
229 xind_l = xIndexFromCj<NbnxnLayout::Simd4xN>(cjFromCi<NbnxnLayout::Simd4xN>(gridj->cell0) + jclusterLast);
231 jx_S = load<SimdReal>(x_j +xind_l + 0*STRIDE_S);
232 jy_S = load<SimdReal>(x_j +xind_l + 1*STRIDE_S);
233 jz_S = load<SimdReal>(x_j +xind_l + 2*STRIDE_S);
235 /* Calculate distance */
236 dx_S0 = load<SimdReal>(x_ci_simd + 0*GMX_SIMD_REAL_WIDTH) - jx_S;
237 dy_S0 = load<SimdReal>(x_ci_simd + 1*GMX_SIMD_REAL_WIDTH) - jy_S;
238 dz_S0 = load<SimdReal>(x_ci_simd + 2*GMX_SIMD_REAL_WIDTH) - jz_S;
239 dx_S1 = load<SimdReal>(x_ci_simd + 3*GMX_SIMD_REAL_WIDTH) - jx_S;
240 dy_S1 = load<SimdReal>(x_ci_simd + 4*GMX_SIMD_REAL_WIDTH) - jy_S;
241 dz_S1 = load<SimdReal>(x_ci_simd + 5*GMX_SIMD_REAL_WIDTH) - jz_S;
242 dx_S2 = load<SimdReal>(x_ci_simd + 6*GMX_SIMD_REAL_WIDTH) - jx_S;
243 dy_S2 = load<SimdReal>(x_ci_simd + 7*GMX_SIMD_REAL_WIDTH) - jy_S;
244 dz_S2 = load<SimdReal>(x_ci_simd + 8*GMX_SIMD_REAL_WIDTH) - jz_S;
245 dx_S3 = load<SimdReal>(x_ci_simd + 9*GMX_SIMD_REAL_WIDTH) - jx_S;
246 dy_S3 = load<SimdReal>(x_ci_simd + 10*GMX_SIMD_REAL_WIDTH) - jy_S;
247 dz_S3 = load<SimdReal>(x_ci_simd + 11*GMX_SIMD_REAL_WIDTH) - jz_S;
249 /* rsq = dx*dx+dy*dy+dz*dz */
250 rsq_S0 = norm2(dx_S0, dy_S0, dz_S0);
251 rsq_S1 = norm2(dx_S1, dy_S1, dz_S1);
252 rsq_S2 = norm2(dx_S2, dy_S2, dz_S2);
253 rsq_S3 = norm2(dx_S3, dy_S3, dz_S3);
255 wco_S0 = (rsq_S0 < rc2_S);
256 wco_S1 = (rsq_S1 < rc2_S);
257 wco_S2 = (rsq_S2 < rc2_S);
258 wco_S3 = (rsq_S3 < rc2_S);
260 wco_any_S01 = wco_S0 || wco_S1;
261 wco_any_S23 = wco_S2 || wco_S3;
262 wco_any_S = wco_any_S01 || wco_any_S23;
264 InRange = anyTrue(wco_any_S);
266 *numDistanceChecks += 4*GMX_SIMD_REAL_WIDTH;
274 if (jclusterFirst <= jclusterLast)
276 for (int jcluster = jclusterFirst; jcluster <= jclusterLast; jcluster++)
278 /* Store cj and the interaction mask */
279 nbl->cj[nbl->ncj].cj = cjFromCi<NbnxnLayout::Simd4xN>(gridj->cell0) + jcluster;
280 nbl->cj[nbl->ncj].excl = get_imask_simd_4xn(excludeSubDiagonal, icluster, jcluster);
283 /* Increase the closing index in i super-cell list */
284 nbl->ci[nbl->nci].cj_ind_end = nbl->ncj;