Expand documentation for nbnxn_put_on_grid()
authorBerk Hess <hess@kth.se>
Wed, 23 Oct 2019 08:27:54 +0000 (10:27 +0200)
committerBerk Hess <hess@kth.se>
Wed, 23 Oct 2019 08:30:36 +0000 (10:30 +0200)
Change-Id: Id6349848d9d03215632cf23e59a556f4e40327fb

src/gromacs/nbnxm/gridset.cpp
src/gromacs/nbnxm/nbnxm.cpp
src/gromacs/nbnxm/nbnxm.h

index b4d76234bf26f5cf65506dcc1190c8194588c8d9..4f476b4e02e283246ecea2cc55e48dd823cbffc5 100644 (file)
@@ -133,7 +133,6 @@ void GridSet::setLocalAtomOrder()
     }
 }
 
-// TODO: Move to gridset.cpp
 void GridSet::putOnGrid(const matrix                    box,
                         const int                       gridIndex,
                         const rvec                      lowerCorner,
index 52d2fecb1305d7c63aabb81608c7416756aa664f..29b68c28b8f1f40fd791799e50411ce76d301938 100644 (file)
@@ -56,7 +56,7 @@
 
 void nbnxn_put_on_grid(nonbonded_verlet_t             *nb_verlet,
                        const matrix                    box,
-                       int                             ddZone,
+                       int                             gridIndex,
                        const rvec                      lowerCorner,
                        const rvec                      upperCorner,
                        const gmx::UpdateGroupsCog     *updateGroupsCog,
@@ -68,7 +68,7 @@ void nbnxn_put_on_grid(nonbonded_verlet_t             *nb_verlet,
                        int                             numAtomsMoved,
                        const int                      *move)
 {
-    nb_verlet->pairSearch_->putOnGrid(box, ddZone, lowerCorner, upperCorner,
+    nb_verlet->pairSearch_->putOnGrid(box, gridIndex, lowerCorner, upperCorner,
                                       updateGroupsCog, atomStart, atomEnd, atomDensity,
                                       atomInfo, x, numAtomsMoved, move,
                                       nb_verlet->nbat.get());
index 9e58003db287e4a7cfbc4ed3a521ce01b9491dad..c93ea94b00f941fe58e9a7b5398328a76cabe30d 100644 (file)
@@ -449,16 +449,34 @@ init_nb_verlet(const gmx::MDLogger     &mdlog,
 /*! \brief Put the atoms on the pair search grid.
  *
  * Only atoms atomStart to atomEnd in x are put on the grid.
- * The atom_density is used to determine the grid size.
- * When atomDensity<=0, the density is determined from atomEnd-atomStart and the corners.
- * With domain decomposition part of the n particles might have migrated,
- * but have not been removed yet. This count is given by nmoved.
- * When move[i] < 0 particle i has migrated and will not be put on the grid.
- * Without domain decomposition move will be NULL.
+ * When \p updateGroupsCog != nullptr, atoms are put on the grid
+ * based on the center of geometry of the group they belong to.
+ * Atoms or COGs of groups should be within the bounding box provided,
+ * this is checked in debug builds when not using update groups.
+ * The atom density is used to determine the grid size when \p gridIndex = 0.
+ * When \p atomDensity <= 0, the density is determined from atomEnd-atomStart
+ * and the bounding box corners.
+ * With domain decomposition, part of the atoms might have migrated,
+ * but have not been removed yet. This count is given by \p numAtomsMoved.
+ * When \p move[i] < 0 particle i has migrated and will not be put on the grid.
+ *
+ * \param[in,out] nb_verlet    The non-bonded object
+ * \param[in]     box          Box used for periodic distance calculations
+ * \param[in]     gridIndex    The index of the grid to spread to, always 0 except with test particle insertion
+ * \param[in]     lowerCorner  Atom groups to be gridded should have coordinates >= this corner
+ * \param[in]     upperCorner  Atom groups to be gridded should have coordinates <= this corner
+ * \param[in]     updateGroupsCog  Centers of geometry for update groups, pass nullptr when not using update groups
+ * \param[in]     atomStart    Start of atom range to grid
+ * \param[in]     atomEnd      End of atom range to grid
+ * \param[in]     atomDensity  An estimate of the atom density, used for peformance optimization and only with \p gridIndex = 0
+ * \param[in]     atomInfo     Atom information flags
+ * \param[in]     x            Coordinates for atoms to grid
+ * \param[in]     numAtomsMoved  The number of atoms that will move to another domain, pass 0 without DD
+ * \param[in]     move         Move flags for atoms, pass nullptr without DD
  */
 void nbnxn_put_on_grid(nonbonded_verlet_t             *nb_verlet,
                        const matrix                    box,
-                       int                             ddZone,
+                       int                             gridIndex,
                        const rvec                      lowerCorner,
                        const rvec                      upperCorner,
                        const gmx::UpdateGroupsCog     *updateGroupsCog,