Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2014-2023 The plumed team
3 : (see the PEOPLE file at the root of the distribution for a list of names)
4 :
5 : See http://www.plumed.org for more information.
6 :
7 : This file is part of plumed, version 2.
8 :
9 : plumed is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : plumed 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
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with plumed. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 : #ifndef __PLUMED_tools_LinkCells_h
23 : #define __PLUMED_tools_LinkCells_h
24 :
25 : #include <limits>
26 : #include <vector>
27 : #include <array>
28 : #include "View.h"
29 : #include "Vector.h"
30 : #include "Pbc.h"
31 :
32 : namespace PLMD {
33 :
34 : class Communicator;
35 :
36 : /// \ingroup TOOLBOX
37 : /// A class for doing link cells
38 82 : class LinkCells {
39 : public:
40 : ///This struct contains the configuration of cells
41 : struct CellCollection {
42 : /// The start of each block corresponding to each link cell
43 : std::vector<unsigned> lcell_starts;
44 : /// The number of atoms in each link cell
45 : std::vector<unsigned> lcell_tots;
46 : /// The atoms ordered by link cells
47 : std::vector<unsigned> lcell_lists;
48 : ///return the sum of the number of atoms in the n cells with most atoms
49 : unsigned getMaximimumCombination(unsigned numCells=27) const;
50 : ///returns a const view of the on indexes the given cell
51 16003445 : inline View<const unsigned> getCellIndexes(const unsigned cellno) const {
52 16003445 : return {lcell_lists.data()+lcell_starts[cellno],lcell_tots[cellno]};
53 : // this will make the debug test fail in the case of
54 : // a series of empty cells in the end of the array,
55 : // because lcell_starts[cellno] will be==cell_lists.size()
56 : // the option abuve have the same problem, but it is masked
57 : // with pointer arithmetic
58 : // and should be ignored by the size passed being 0
59 : //return {&lcell_lists[lcell_starts[cellno]],lcell_tots[cellno]};
60 : }
61 : };
62 : private:
63 : /// Minimum number of atoms per thread when finding the cells of the atoms
64 : static constexpr unsigned minAtomsPerThread=4096;
65 : /// Maximum total number of cells
66 : static constexpr double maxNumberOfCells=1048576;
67 : /// Check that the link cells were set up correctly
68 : bool cutoffwasset=false;
69 : /// Are there periodic boundary conditions setup
70 : bool nopbc=false;
71 : /// The cutoff to use for the sizes of the cells
72 : double link_cutoff=0.0;
73 : /// The pbc we are using for link cells
74 : Pbc mypbc;
75 : /// The location of the origin if we are not using periodic boundary conditions
76 : Vector origin;
77 : /// The number of cells in each direction
78 : std::array<unsigned,3> ncells{0,0,0};
79 : /// The number of cells to stride through to get the link cells
80 : std::array<unsigned,3> nstride{1,0,0};
81 : /// Work vector with the list of cells each atom is inside
82 : std::vector<unsigned> allcells;
83 : ///The collection of indexes per cell created by buildCellLists
84 : CellCollection innerCollection;
85 : void createCells(const PLMD::Tensor&);
86 : /// Find the index of the cell of a position, without checks:
87 : /// returns getNumberOfCells() if the position is not in a cell
88 : unsigned findCellUnchecked( Vector pos ) const noexcept;
89 : public:
90 : /// The communicator is not used: the cells are built locally on every process
91 : /// (the argument is kept for backward compatibility)
92 : explicit LinkCells( Communicator& comm );
93 : /// Have the link cells been enabled
94 : bool enabled() const ;
95 : /// Set the value of the cutoff
96 : void setCutoff( double lcut );
97 : /// Get the value of the cutoff
98 : double getCutoff() const ;
99 : /// Get the total number of link cells
100 : unsigned getNumberOfCells() const ;
101 : ///Get the number of cells per dimension
102 : const std::array<unsigned,3>& getCellLimits() const ;
103 : /// Get the number of atoms in the cell that contains the most atoms
104 : unsigned getMaxInCell() const ;
105 : /// Sets up the cells ignoring the pbcs
106 : ///
107 : /// This creates an orthogonal box that encloses all the atoms
108 : void setupCells( View<const Vector> pos);
109 : template <typename T>
110 : inline void setupCells( const T& pos ) {
111 : setupCells(make_const_view(pos));
112 : }
113 : /// Sets up the cells using the box stored in the Pbc object
114 : void setupCells( const Pbc& pbc );
115 : ///setups the cells
116 : ///
117 : /// Usually it only uses the pbcs
118 : /// if the pbc box is null (three zero vectors) it builds an orthorhombic box from the given positions
119 : void setupCells( View<const Vector> pos,
120 : const Pbc& pbc );
121 : template <typename T>
122 : inline void setupCells( const T& pos,
123 : const Pbc& pbc ) {
124 : setupCells(make_const_view(pos),pbc);
125 : }
126 : /// Build the link cell lists
127 : void buildCellLists( View<const Vector> pos,
128 : View<const unsigned> indices,
129 : const Pbc& pbc );
130 : template <typename T,typename Y>
131 : inline void buildCellLists( const T& pos,
132 : const Y& indices,
133 : const Pbc& pbc ) {
134 : buildCellLists(make_const_view(pos), make_const_view(indices),pbc);
135 : }
136 : // Reset the passed collection on the positions and the passed indexes, needs setupCells to be launched in advance
137 : void resetCollection( CellCollection& collection,
138 : View<const Vector> pos,
139 : View<const unsigned> indices );
140 : // Returns a collection on the positions and the passed indexes, needs setupCells to be launched in advance
141 : CellCollection getCollection(View<const Vector> pos,
142 : View<const unsigned> indices );
143 : template <typename T,typename Y>
144 : inline CellCollection getCollection(const T& pos,
145 : const Y& indices) {
146 : return getCollection(make_const_view(pos),
147 : make_const_view(indices));
148 : }
149 : /// Take three indices and return the index of the corresponding cell
150 : unsigned convertIndicesToIndex( std::array<unsigned,3> cellCoord) const;
151 : /// Find the cell index in which this position is contained
152 : unsigned findCell( const Vector& pos ) const ;
153 : /// Find the cell in which this position is contained
154 : std::array<unsigned,3> findMyCell( Vector pos ) const ;
155 : /// Returns the "coordinate" of the cell index, assumes that cellIndex < getNumberOfCells()
156 : std::array<unsigned,3> findMyCell( unsigned cellIndex ) const ;
157 : /// Get the list of cells we need to surround the a particular cell
158 : ///
159 : /// usePbc toggles the use of pbcs for the input cell (celln)
160 : void addRequiredCells( const std::array<unsigned,3>& celn,
161 : unsigned& ncells_required,
162 : std::vector<unsigned>& cells_required,
163 : bool usePbc=true) const ;
164 : /// Retrieve the atoms in a list of cells
165 : void retrieveAtomsInCells( unsigned ncells_required,
166 : View<const unsigned> cells_required,
167 : unsigned& natomsper,
168 : std::vector<unsigned>& atoms,
169 : unsigned avoidIndex=std::numeric_limits<unsigned>::max()) const ;
170 : /// Retrieve the atoms we need to consider
171 : void retrieveNeighboringAtoms( const Vector& pos,
172 : std::vector<unsigned>& cell_list,
173 : unsigned& natomsper,
174 : std::vector<unsigned>& atoms ) const ;
175 : /// Create a neighbour list for the specified input atoms
176 : ///
177 : /// @param nat number of atoms
178 : /// @param pos position of the "central atoms"
179 : /// @param ind indexes of the central atoms
180 : /// @param tind temporary indexes
181 : /// @param neigh_pos position of the candidate neighbors
182 : /// @param neigh_ind indexes of the candidate neighbors
183 : /// @param pbc
184 : /// @param natoms_per_list *output* number of atom per list (mapped on tind) `nlist[tind[i]]=natoms[i]`
185 : /// @param nlist *output* list of the central atoms plus their neighbors
186 : /// @param applyCutoff if true only the candidates whose (minimum image) distance from the
187 : /// central atom is not larger than the link cell cutoff are kept. The cells are at
188 : /// least as wide as the cutoff, so without this the lists hold all the atoms of the
189 : /// 27 cells around the central atom, up to ~2-3 cutoffs away: with a short cutoff
190 : /// in a dense system most of them are beyond it. Use it only when the quantity
191 : /// computed from the list is exactly zero beyond the cutoff.
192 : void createNeighborList( unsigned nat,
193 : View<const Vector> pos,
194 : View<const unsigned> ind,
195 : View<const unsigned> tind,
196 : View<const Vector> neigh_pos,
197 : View<const unsigned> neigh_ind,
198 : const Pbc& pbc,
199 : unsigned& natoms_per_list,
200 : std::vector<std::size_t>& nlist,
201 : const bool applyCutoff=false ) ;
202 : };
203 :
204 : inline
205 : bool LinkCells::enabled() const {
206 : return cutoffwasset;
207 : }
208 :
209 : inline
210 : unsigned LinkCells::getNumberOfCells() const {
211 15953118 : return ncells[0]*ncells[1]*ncells[2];
212 : }
213 :
214 : inline
215 : const std::array<unsigned,3>& LinkCells::getCellLimits() const {
216 : return ncells;
217 : }
218 : }
219 :
220 : #endif
|
