Enable TPI with the Verlet cut-off scheme
[gromacs.git] / src / gromacs / nbnxm / gridset.h
blob142a1f28b4126b53a18e6638c59d5dc197369734
1 /*
2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 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.
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 /*! \internal \file
38 * \brief
39 * Declares the GridSet class.
41 * This class holds the grids for the local and non-local domain decomposition
42 * zones, as well as the cell and atom data that covers all grids.
44 * \author Berk Hess <hess@kth.se>
45 * \ingroup module_nbnxm
48 #ifndef GMX_NBNXM_GRIDSET_H
49 #define GMX_NBNXM_GRIDSET_H
51 #include <memory>
52 #include <vector>
54 #include "gromacs/math/vec.h"
55 #include "gromacs/math/vectypes.h"
56 #include "gromacs/utility/alignedallocator.h"
57 #include "gromacs/utility/arrayref.h"
59 #include "grid.h"
60 #include "gridsetdata.h"
63 struct nbnxn_atomdata_t;
64 enum class PairlistType;
66 namespace gmx
68 class UpdateGroupsCog;
71 namespace Nbnxm
74 /*! \internal
75 * \brief Holds a set of search grids for the local + non-local DD zones
77 * The are three different possible setups:
78 * - a single grid, this is the standard case without domain decomposition
79 * - one grid for each domain decomposition zone
80 * - with test particle insertion there are two grids, one for the system
81 * to insert in and one for the molecule that is inserted
83 class GridSet
85 public:
86 /*! \internal
87 * \brief Description of the domain setup: PBC and the connections between domains
89 struct DomainSetup
91 //! Constructor, without DD \p numDDCells and \p ddZones should be nullptr
92 DomainSetup(int ePBC,
93 bool doTestParticleInsertion,
94 const ivec *numDDCells,
95 const gmx_domdec_zones_t *ddZones);
97 //! The type of PBC
98 int ePBC;
99 //! Tells whether we are doing test-particle insertion
100 bool doTestParticleInsertion;
101 //! Are there multiple domains?
102 bool haveMultipleDomains;
103 //! Are there multiple domains along each dimension?
104 std::array<bool, DIM> haveMultipleDomainsPerDim;
105 //! The domain decomposition zone setup
106 const gmx_domdec_zones_t *zones;
109 //! Constructs a grid set for 1 or multiple DD zones, when numDDCells!=nullptr
110 GridSet(int ePBC,
111 bool doTestParticleInsertion,
112 const ivec *numDDCells,
113 const gmx_domdec_zones_t *ddZones,
114 PairlistType pairlistType,
115 bool haveFep,
116 int numThreads,
117 gmx::PinningPolicy pinningPolicy);
119 //! Puts the atoms on the grid with index \p gridIndex and copies the coordinates to \p nbat
120 void putOnGrid(const matrix box,
121 int gridIndex,
122 const rvec lowerCorner,
123 const rvec upperCorner,
124 const gmx::UpdateGroupsCog *updateGroupsCog,
125 int atomStart,
126 int atomEnd,
127 real atomDensity,
128 gmx::ArrayRef<const int> atomInfo,
129 gmx::ArrayRef<const gmx::RVec> x,
130 int numAtomsMoved,
131 const int *move,
132 nbnxn_atomdata_t *nbat);
134 //! Returns the domain setup
135 DomainSetup domainSetup() const
137 return domainSetup_;
140 //! Returns the total number of atoms in the grid set, including padding
141 int numGridAtomsTotal() const
143 return grids_.back().atomIndexEnd();
146 //! Returns the number of local real atoms, i.e. without padded atoms
147 int numRealAtomsLocal() const
149 return numRealAtomsLocal_;
152 //! Returns the number of total real atoms, i.e. without padded atoms
153 int numRealAtomsTotal() const
155 return numRealAtomsTotal_;
158 //! Returns the atom order on the grid for the local atoms
159 gmx::ArrayRef<const int> getLocalAtomorder() const
161 /* Return the atom order for the home cell (index 0) */
162 const int numIndices = grids_[0].atomIndexEnd() - grids_[0].firstAtomInColumn(0);
164 return gmx::constArrayRefFromArray(atomIndices().data(), numIndices);
167 //! Sets the order of the local atoms to the order grid atom ordering
168 void setLocalAtomOrder();
170 //! Returns the list of grids
171 gmx::ArrayRef<const Grid> grids() const
173 return grids_;
176 //! Returns the grid atom indices covering all grids
177 gmx::ArrayRef<const int> cells() const
179 return gridSetData_.cells;
182 //! Returns the grid atom indices covering all grids
183 gmx::ArrayRef<const int> atomIndices() const
185 return gridSetData_.atomIndices;
188 //! Returns whether we have perturbed non-bonded interactions
189 bool haveFep() const
191 return haveFep_;
194 //! Returns the unit cell in \p box
195 void getBox(matrix box) const
197 copy_mat(box_, box);
200 //! Returns the maximum number of columns across all grids
201 int numColumnsMax() const
203 return numColumnsMax_;
206 //! Sets the maximum number of columns across all grids
207 void setNumColumnsMax(int numColumnsMax)
209 numColumnsMax_ = numColumnsMax;
212 private:
213 /* Data members */
214 //! The domain setup
215 DomainSetup domainSetup_;
216 //! The search grids
217 std::vector<Grid> grids_;
218 //! The cell and atom index data which runs over all grids
219 GridSetData gridSetData_;
220 //! Tells whether we have perturbed non-bonded interactions
221 bool haveFep_;
222 //! The periodic unit-cell
223 matrix box_;
224 //! The number of local real atoms, i.e. without padded atoms, local atoms: 0 to numAtomsLocal_
225 int numRealAtomsLocal_;
226 //! The total number of real atoms, i.e. without padded atoms
227 int numRealAtomsTotal_;
228 //! Working data for constructing a single grid, one entry per thread
229 std::vector<GridWork> gridWork_;
230 //! Maximum number of columns across all grids
231 int numColumnsMax_;
235 } // namespace Nbnxm
237 #endif