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 >= 2*NBNXN_CPU_CLUSTER_I_SIZE
37 #define STRIDE_S (GMX_SIMD_REAL_WIDTH/2)
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_2xnn(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::Simd2xNN
>(ci
);
54 store(x_ci_simd
+ 0*GMX_SIMD_REAL_WIDTH
, loadU1DualHsimd(x
+ ia
+ 0*STRIDE_S
+ 0) + SimdReal(shx
) );
55 store(x_ci_simd
+ 1*GMX_SIMD_REAL_WIDTH
, loadU1DualHsimd(x
+ ia
+ 1*STRIDE_S
+ 0) + SimdReal(shy
) );
56 store(x_ci_simd
+ 2*GMX_SIMD_REAL_WIDTH
, loadU1DualHsimd(x
+ ia
+ 2*STRIDE_S
+ 0) + SimdReal(shz
) );
57 store(x_ci_simd
+ 3*GMX_SIMD_REAL_WIDTH
, loadU1DualHsimd(x
+ ia
+ 0*STRIDE_S
+ 2) + SimdReal(shx
) );
58 store(x_ci_simd
+ 4*GMX_SIMD_REAL_WIDTH
, loadU1DualHsimd(x
+ ia
+ 1*STRIDE_S
+ 2) + SimdReal(shy
) );
59 store(x_ci_simd
+ 5*GMX_SIMD_REAL_WIDTH
, loadU1DualHsimd(x
+ ia
+ 2*STRIDE_S
+ 2) + SimdReal(shz
) );
62 /* SIMD code for checking and adding cluster-pairs to the list using coordinates in packed format.
64 * Checks bouding box distances and possibly atom pair distances.
65 * This is an accelerated version of make_cluster_list_simple.
67 * \param[in] gridj The j-grid
68 * \param[in,out] nbl The pair-list to store the cluster pairs in
69 * \param[in] icluster The index of the i-cluster
70 * \param[in] firstCell The first cluster in the j-range, using i-cluster size indexing
71 * \param[in] lastCell The last cluster in the j-range, using i-cluster size indexing
72 * \param[in] excludeSubDiagonal Exclude atom pairs with i-index > j-index
73 * \param[in] x_j Coordinates for the j-atom, in SIMD packed format
74 * \param[in] rlist2 The squared list cut-off
75 * \param[in] rbb2 The squared cut-off for putting cluster-pairs in the list based on bounding box distance only
76 * \param[in,out] numDistanceChecks The number of distance checks performed
79 makeClusterListSimd2xnn(const nbnxn_grid_t
* gridj
,
80 nbnxn_pairlist_t
* nbl
,
84 bool excludeSubDiagonal
,
85 const real
* gmx_restrict x_j
,
88 int * gmx_restrict numDistanceChecks
)
91 const real
* gmx_restrict x_ci_simd
= nbl
->work
->x_ci_simd
;
92 const nbnxn_bb_t
* gmx_restrict bb_ci
= nbl
->work
->bb_ci
;
94 SimdReal jx_S
, jy_S
, jz_S
;
96 SimdReal dx_S0
, dy_S0
, dz_S0
;
97 SimdReal dx_S2
, dy_S2
, dz_S2
;
112 int jclusterFirst
= cjFromCi
<NbnxnLayout::Simd2xNN
>(firstCell
);
113 int jclusterLast
= cjFromCi
<NbnxnLayout::Simd2xNN
>(lastCell
);
114 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");
116 rc2_S
= SimdReal(rlist2
);
119 while (!InRange
&& jclusterFirst
<= jclusterLast
)
121 #if NBNXN_SEARCH_BB_SIMD4
122 d2
= subc_bb_dist2_simd4(0, bb_ci
, jclusterFirst
, gridj
->bbj
);
124 d2
= subc_bb_dist2(0, bb_ci
, jclusterFirst
, gridj
->bbj
);
126 *numDistanceChecks
+= 2;
128 /* Check if the distance is within the distance where
129 * we use only the bounding box distance rbb,
130 * or within the cut-off and there is at least one atom pair
131 * within the cut-off.
137 else if (d2
< rlist2
)
139 xind_f
= xIndexFromCj
<NbnxnLayout::Simd2xNN
>(cjFromCi
<NbnxnLayout::Simd2xNN
>(gridj
->cell0
) + jclusterFirst
);
141 jx_S
= loadDuplicateHsimd(x_j
+ xind_f
+ 0*STRIDE_S
);
142 jy_S
= loadDuplicateHsimd(x_j
+ xind_f
+ 1*STRIDE_S
);
143 jz_S
= loadDuplicateHsimd(x_j
+ xind_f
+ 2*STRIDE_S
);
145 /* Calculate distance */
146 dx_S0
= load
<SimdReal
>(x_ci_simd
+ 0*GMX_SIMD_REAL_WIDTH
) - jx_S
;
147 dy_S0
= load
<SimdReal
>(x_ci_simd
+ 1*GMX_SIMD_REAL_WIDTH
) - jy_S
;
148 dz_S0
= load
<SimdReal
>(x_ci_simd
+ 2*GMX_SIMD_REAL_WIDTH
) - jz_S
;
149 dx_S2
= load
<SimdReal
>(x_ci_simd
+ 3*GMX_SIMD_REAL_WIDTH
) - jx_S
;
150 dy_S2
= load
<SimdReal
>(x_ci_simd
+ 4*GMX_SIMD_REAL_WIDTH
) - jy_S
;
151 dz_S2
= load
<SimdReal
>(x_ci_simd
+ 5*GMX_SIMD_REAL_WIDTH
) - jz_S
;
153 /* rsq = dx*dx+dy*dy+dz*dz */
154 rsq_S0
= norm2(dx_S0
, dy_S0
, dz_S0
);
155 rsq_S2
= norm2(dx_S2
, dy_S2
, dz_S2
);
157 wco_S0
= (rsq_S0
< rc2_S
);
158 wco_S2
= (rsq_S2
< rc2_S
);
160 wco_any_S
= wco_S0
|| wco_S2
;
162 InRange
= anyTrue(wco_any_S
);
164 *numDistanceChecks
+= 2*GMX_SIMD_REAL_WIDTH
;
177 while (!InRange
&& jclusterLast
> jclusterFirst
)
179 #if NBNXN_SEARCH_BB_SIMD4
180 d2
= subc_bb_dist2_simd4(0, bb_ci
, jclusterLast
, gridj
->bbj
);
182 d2
= subc_bb_dist2(0, bb_ci
, jclusterLast
, gridj
->bbj
);
184 *numDistanceChecks
+= 2;
186 /* Check if the distance is within the distance where
187 * we use only the bounding box distance rbb,
188 * or within the cut-off and there is at least one atom pair
189 * within the cut-off.
195 else if (d2
< rlist2
)
197 xind_l
= xIndexFromCj
<NbnxnLayout::Simd2xNN
>(cjFromCi
<NbnxnLayout::Simd2xNN
>(gridj
->cell0
) + jclusterLast
);
199 jx_S
= loadDuplicateHsimd(x_j
+ xind_l
+ 0*STRIDE_S
);
200 jy_S
= loadDuplicateHsimd(x_j
+ xind_l
+ 1*STRIDE_S
);
201 jz_S
= loadDuplicateHsimd(x_j
+ xind_l
+ 2*STRIDE_S
);
203 /* Calculate distance */
204 dx_S0
= load
<SimdReal
>(x_ci_simd
+ 0*GMX_SIMD_REAL_WIDTH
) - jx_S
;
205 dy_S0
= load
<SimdReal
>(x_ci_simd
+ 1*GMX_SIMD_REAL_WIDTH
) - jy_S
;
206 dz_S0
= load
<SimdReal
>(x_ci_simd
+ 2*GMX_SIMD_REAL_WIDTH
) - jz_S
;
207 dx_S2
= load
<SimdReal
>(x_ci_simd
+ 3*GMX_SIMD_REAL_WIDTH
) - jx_S
;
208 dy_S2
= load
<SimdReal
>(x_ci_simd
+ 4*GMX_SIMD_REAL_WIDTH
) - jy_S
;
209 dz_S2
= load
<SimdReal
>(x_ci_simd
+ 5*GMX_SIMD_REAL_WIDTH
) - jz_S
;
211 /* rsq = dx*dx+dy*dy+dz*dz */
212 rsq_S0
= norm2(dx_S0
, dy_S0
, dz_S0
);
213 rsq_S2
= norm2(dx_S2
, dy_S2
, dz_S2
);
215 wco_S0
= (rsq_S0
< rc2_S
);
216 wco_S2
= (rsq_S2
< rc2_S
);
218 wco_any_S
= wco_S0
|| wco_S2
;
220 InRange
= anyTrue(wco_any_S
);
222 *numDistanceChecks
+= 2*GMX_SIMD_REAL_WIDTH
;
230 if (jclusterFirst
<= jclusterLast
)
232 for (int jcluster
= jclusterFirst
; jcluster
<= jclusterLast
; jcluster
++)
234 /* Store cj and the interaction mask */
235 nbl
->cj
[nbl
->ncj
].cj
= cjFromCi
<NbnxnLayout::Simd2xNN
>(gridj
->cell0
) + jcluster
;
236 nbl
->cj
[nbl
->ncj
].excl
= get_imask_simd_2xnn(excludeSubDiagonal
, icluster
, jcluster
);
239 /* Increase the closing index in i super-cell list */
240 nbl
->ci
[nbl
->nci
].cj_ind_end
= nbl
->ncj
;