src/gromacs/nbnxm/pairlist_simd_2xmm.h

   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 /*! \internal \file
  37  *
  38  * \brief
  39  * Declares inline-friendly code for making 2xNN pairlists
  40  *
  41  * \author Berk Hess <hess@kth.se>
  42  * \ingroup module_nbnxm
  43  */
  44
  45
  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;
  50
  51 //! Copies PBC shifted i-cell packed atom coordinates to working array
  52 static inline void icell_set_x_simd_2xnn(int  ci,
  53                                          real shx,
  54                                          real shy,
  55                                          real shz,
  56                                          int gmx_unused        stride,
  57                                          const real*           x,
  58                                          NbnxnPairlistCpuWork* work)
  59 {
  60     int   ia;
  61     real* x_ci_simd = work->iClusterData.xSimd.data();
  62
  63     ia = xIndexFromCi<NbnxnLayout::Simd2xNN>(ci);
  64
  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));
  77 }
  78
  79 /*! \brief SIMD code for checking and adding cluster-pairs to the list using coordinates in packed format.
  80  *
  81  * Checks bouding box distances and possibly atom pair distances.
  82  * This is an accelerated version of make_cluster_list_simple.
  83  *
  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
  94  */
  95 static inline void makeClusterListSimd2xnn(const Grid&       jGrid,
  96                                            NbnxnPairlistCpu* nbl,
  97                                            int               icluster,
  98                                            int               firstCell,
  99                                            int               lastCell,
 100                                            bool              excludeSubDiagonal,
 101                                            const real* gmx_restrict x_j,
 102                                            real                     rlist2,
 103                                            float                    rbb2,
 104                                            int* gmx_restrict numDistanceChecks)
 105 {
 106     using namespace gmx;
 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();
 109
 110     SimdReal jx_S, jy_S, jz_S;
 111
 112     SimdReal dx_S0, dy_S0, dz_S0;
 113     SimdReal dx_S2, dy_S2, dz_S2;
 114
 115     SimdReal rsq_S0;
 116     SimdReal rsq_S2;
 117
 118     SimdBool wco_S0;
 119     SimdBool wco_S2;
 120     SimdBool wco_any_S;
 121
 122     SimdReal rc2_S;
 123
 124     gmx_bool InRange;
 125     float    d2;
 126     int      xind_f, xind_l;
 127
 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");
 133
 134     rc2_S = SimdReal(rlist2);
 135
 136     InRange = FALSE;
 137     while (!InRange && jclusterFirst <= jclusterLast)
 138     {
 139         d2 = clusterBoundingBoxDistance2(bb_ci[0], jGrid.jBoundingBoxes()[jclusterFirst]);
 140         *numDistanceChecks += 2;
 141
 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.
 146          */
 147         if (d2 < rbb2)
 148         {
 149             InRange = TRUE;
 150         }
 151         else if (d2 < rlist2)
 152         {
 153             xind_f = xIndexFromCj<NbnxnLayout::Simd2xNN>(
 154                     cjFromCi<NbnxnLayout::Simd2xNN, 0>(jGrid.cellOffset()) + jclusterFirst);
 155
 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);
 159
 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;
 167
 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);
 171
 172             wco_S0 = (rsq_S0 < rc2_S);
 173             wco_S2 = (rsq_S2 < rc2_S);
 174
 175             wco_any_S = wco_S0 || wco_S2;
 176
 177             InRange = anyTrue(wco_any_S);
 178
 179             *numDistanceChecks += 2 * GMX_SIMD_REAL_WIDTH;
 180         }
 181         if (!InRange)
 182         {
 183             jclusterFirst++;
 184         }
 185     }
 186     if (!InRange)
 187     {
 188         return;
 189     }
 190
 191     InRange = FALSE;
 192     while (!InRange && jclusterLast > jclusterFirst)
 193     {
 194         d2 = clusterBoundingBoxDistance2(bb_ci[0], jGrid.jBoundingBoxes()[jclusterLast]);
 195         *numDistanceChecks += 2;
 196
 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.
 201          */
 202         if (d2 < rbb2)
 203         {
 204             InRange = TRUE;
 205         }
 206         else if (d2 < rlist2)
 207         {
 208             xind_l = xIndexFromCj<NbnxnLayout::Simd2xNN>(
 209                     cjFromCi<NbnxnLayout::Simd2xNN, 0>(jGrid.cellOffset()) + jclusterLast);
 210
 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);
 214
 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;
 222
 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);
 226
 227             wco_S0 = (rsq_S0 < rc2_S);
 228             wco_S2 = (rsq_S2 < rc2_S);
 229
 230             wco_any_S = wco_S0 || wco_S2;
 231
 232             InRange = anyTrue(wco_any_S);
 233
 234             *numDistanceChecks += 2 * GMX_SIMD_REAL_WIDTH;
 235         }
 236         if (!InRange)
 237         {
 238             jclusterLast--;
 239         }
 240     }
 241
 242     if (jclusterFirst <= jclusterLast)
 243     {
 244         for (int jcluster = jclusterFirst; jcluster <= jclusterLast; jcluster++)
 245         {
 246             /* Store cj and the interaction mask */
 247             nbnxn_cj_t cjEntry;
 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);
 251         }
 252         /* Increase the closing index in the i list */
 253         nbl->ci.back().cj_ind_end = nbl->cj.size();
 254     }
 255 }