Make nbnxm a proper module
[alexxy/gromacs.git] / src / gromacs / nbnxm / kernels_simd_2xmm / kernel_prune.cpp
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 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.
8  *
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.
13  *
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.
18  *
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.
23  *
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.
31  *
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.
34  */
35
36 #include "gmxpre.h"
37
38 #include "kernel_prune.h"
39
40 #include "gromacs/nbnxm/atomdata.h"
41 #include "gromacs/nbnxm/nbnxm_simd.h"
42 #include "gromacs/nbnxm/pairlist.h"
43 #include "gromacs/utility/gmxassert.h"
44
45 #ifdef GMX_NBNXN_SIMD_2XNN
46 #define GMX_SIMD_J_UNROLL_SIZE 2
47 #include "kernel_common.h"
48 #endif
49
50 /* Prune a single nbnxn_pairtlist_t entry with distance rlistInner */
51 void
52 nbnxn_kernel_prune_2xnn(NbnxnPairlistCpu *         nbl,
53                         const nbnxn_atomdata_t *   nbat,
54                         const rvec * gmx_restrict  shift_vec,
55                         real                       rlistInner)
56 {
57 #ifdef GMX_NBNXN_SIMD_2XNN
58     using namespace gmx;
59
60     /* We avoid push_back() for efficiency reasons and resize after filling */
61     nbl->ci.resize(nbl->ciOuter.size());
62     nbl->cj.resize(nbl->cjOuter.size());
63
64     const nbnxn_ci_t * gmx_restrict ciOuter  = nbl->ciOuter.data();
65     nbnxn_ci_t       * gmx_restrict ciInner  = nbl->ci.data();
66
67     const nbnxn_cj_t * gmx_restrict cjOuter  = nbl->cjOuter.data();
68     nbnxn_cj_t       * gmx_restrict cjInner  = nbl->cj.data();
69
70     const real       * gmx_restrict shiftvec = shift_vec[0];
71     const real       * gmx_restrict x        = nbat->x().data();
72
73     const SimdReal                  rlist2_S(rlistInner*rlistInner);
74
75     /* Initialize the new list count as empty and add pairs that are in range */
76     int       nciInner = 0;
77     int       ncjInner = 0;
78     const int nciOuter = nbl->ciOuter.size();
79     for (int i = 0; i < nciOuter; i++)
80     {
81         const nbnxn_ci_t * gmx_restrict ciEntry = &ciOuter[i];
82
83         /* Copy the original list entry to the pruned entry */
84         ciInner[nciInner].ci           = ciEntry->ci;
85         ciInner[nciInner].shift        = ciEntry->shift;
86         ciInner[nciInner].cj_ind_start = ncjInner;
87
88         /* Extract shift data */
89         int      ish     = (ciEntry->shift & NBNXN_CI_SHIFT);
90         int      ish3    = ish*3;
91         int      ci      = ciEntry->ci;
92
93         SimdReal shX_S   = SimdReal(shiftvec[ish3    ]);
94         SimdReal shY_S   = SimdReal(shiftvec[ish3 + 1]);
95         SimdReal shZ_S   = SimdReal(shiftvec[ish3 + 2]);
96
97 #if UNROLLJ <= 4
98         int      scix    = ci*STRIDE*DIM;
99 #else
100         int      scix    = (ci >> 1)*STRIDE*DIM + (ci & 1)*(STRIDE >> 1);
101 #endif
102
103         /* Load i atom data */
104         int      sciy    = scix + STRIDE;
105         int      sciz    = sciy + STRIDE;
106         SimdReal ix_S0   = loadU1DualHsimd(x + scix    ) + shX_S;
107         SimdReal ix_S2   = loadU1DualHsimd(x + scix + 2) + shX_S;
108         SimdReal iy_S0   = loadU1DualHsimd(x + sciy    ) + shY_S;
109         SimdReal iy_S2   = loadU1DualHsimd(x + sciy + 2) + shY_S;
110         SimdReal iz_S0   = loadU1DualHsimd(x + sciz    ) + shZ_S;
111         SimdReal iz_S2   = loadU1DualHsimd(x + sciz + 2) + shZ_S;
112
113         for (int cjind = ciEntry->cj_ind_start; cjind < ciEntry->cj_ind_end; cjind++)
114         {
115             /* j-cluster index */
116             int cj      = cjOuter[cjind].cj;
117
118             /* Atom indices (of the first atom in the cluster) */
119 #if UNROLLJ == STRIDE
120             int aj      = cj*UNROLLJ;
121             int ajx     = aj*DIM;
122 #else
123             int ajx     = (cj >> 1)*DIM*STRIDE + (cj & 1)*UNROLLJ;
124 #endif
125             int ajy     = ajx + STRIDE;
126             int ajz     = ajy + STRIDE;
127
128             /* load j atom coordinates */
129             SimdReal jx_S   = loadDuplicateHsimd(x + ajx);
130             SimdReal jy_S   = loadDuplicateHsimd(x + ajy);
131             SimdReal jz_S   = loadDuplicateHsimd(x + ajz);
132
133             /* Calculate distance */
134             SimdReal dx_S0  = ix_S0 - jx_S;
135             SimdReal dy_S0  = iy_S0 - jy_S;
136             SimdReal dz_S0  = iz_S0 - jz_S;
137             SimdReal dx_S2  = ix_S2 - jx_S;
138             SimdReal dy_S2  = iy_S2 - jy_S;
139             SimdReal dz_S2  = iz_S2 - jz_S;
140
141             /* rsq = dx*dx+dy*dy+dz*dz */
142             SimdReal rsq_S0 = norm2(dx_S0, dy_S0, dz_S0);
143             SimdReal rsq_S2 = norm2(dx_S2, dy_S2, dz_S2);
144
145             /* Do the cut-off check */
146             SimdBool wco_S0 = (rsq_S0 < rlist2_S);
147             SimdBool wco_S2 = (rsq_S2 < rlist2_S);
148
149             wco_S0          = wco_S0 || wco_S2;
150
151             /* Putting the assignment inside the conditional is slower */
152             cjInner[ncjInner] = cjOuter[cjind];
153             if (anyTrue(wco_S0))
154             {
155                 ncjInner++;
156             }
157         }
158
159         if (ncjInner > ciInner[nciInner].cj_ind_start)
160         {
161             ciInner[nciInner].cj_ind_end = ncjInner;
162             nciInner++;
163         }
164     }
165
166     nbl->ci.resize(nciInner);
167     nbl->cj.resize(ncjInner);
168
169 #else  /* GMX_NBNXN_SIMD_2XNN */
170
171     GMX_RELEASE_ASSERT(false, "2xNN kernel called without 2xNN support");
172
173     GMX_UNUSED_VALUE(nbl);
174     GMX_UNUSED_VALUE(nbat);
175     GMX_UNUSED_VALUE(shift_vec);
176     GMX_UNUSED_VALUE(rlistInner);
177
178 #endif /* GMX_NBNXN_SIMD_2XNN */
179 }