1 /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 Copyright (c) 2011-2021 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_core_Atoms_h
23 #define __PLUMED_core_Atoms_h
24
25 #include "tools/Communicator.h"
26 #include "tools/Tensor.h"
27 #include "tools/Units.h"
28 #include "tools/Exception.h"
29 #include "tools/AtomNumber.h"
30 #include "tools/ForwardDecl.h"
31 #include <vector>
32 #include <set>
33 #include <map>
34 #include <string>
35 #include <memory>
36
37 namespace PLMD {
38
39 class MDAtomsBase;
40 class PlumedMain;
41 class ActionAtomistic;
42 class ActionWithVirtualAtom;
43 class Pbc;
44
45 /// Class containing atom related quantities from the MD code.
46 /// IT IS STILL UNDOCUMENTED. IT PROBABLY NEEDS A STRONG CLEANUP
47 class Atoms
48 {
49 friend class ActionAtomistic;
50 friend class ActionWithVirtualAtom;
51 int natoms;
52 std::set<AtomNumber> unique;
53 std::vector<unsigned> uniq_index;
54 /// Map global indexes to local indexes
55 /// E.g. g2l[i] is the position of atom i in the array passed from the MD engine.
56 /// Called "global to local" since originally it was used to map global indexes to local
57 /// ones used in domain decomposition. However, it is now also used for the NAMD-like
58 /// interface, where only a small number of atoms is passed to plumed.
59 std::vector<int> g2l;
60 std::vector<Vector> positions;
61 std::vector<Vector> forces;
62 std::vector<double> masses;
63 std::vector<double> charges;
64 std::vector<ActionWithVirtualAtom*> virtualAtomsActions;
65 Tensor box;
66 ForwardDecl<Pbc> pbc_fwd;
67 Pbc& pbc=*pbc_fwd;
68 Tensor virial;
69 // this is the energy set by each processor:
70 double md_energy;
71 // this is the summed energy:
72 double energy;
73
74 bool dataCanBeSet;
75 bool collectEnergy;
76 bool energyHasBeenSet;
77 unsigned positionsHaveBeenSet;
78 bool massesHaveBeenSet;
79 bool chargesHaveBeenSet;
80 bool boxHasBeenSet;
81 unsigned forcesHaveBeenSet;
82 bool virialHasBeenSet;
83 bool massAndChargeOK;
84 unsigned shuffledAtoms;
85
86 std::map<std::string,std::vector<AtomNumber> > groups;
87
88 void resizeVectors(unsigned);
89
90 std::vector<int> fullList;
91
92 std::unique_ptr<MDAtomsBase> mdatoms;
93
94 PlumedMain & plumed;
95
96 Units MDUnits;
97 Units units;
98
99 bool naturalUnits;
100 bool MDnaturalUnits;
101
102 double timestep;
103 double forceOnEnergy;
104
105 /// if set to true, all the forces in the global array are zeroes
106 /// at every step. It should not be necessary in general, but it is
107 /// for actions accessing to modifyGlobalForce() (e.g. FIT_TO_TEMPLATE).
108 bool zeroallforces;
109
110 double kbT;
111
112 std::vector<ActionAtomistic*> actions;
113 std::vector<int> gatindex;
114
115 bool asyncSent;
116 bool atomsNeeded;
117
118 class DomainDecomposition:
119 public Communicator
120 {
121 public:
122 bool on;
123 bool async;
124
125 std::vector<Communicator::Request> mpi_request_positions;
126 std::vector<Communicator::Request> mpi_request_index;
127
128 std::vector<double> positionsToBeSent;
129 std::vector<double> positionsToBeReceived;
130 std::vector<int> indexToBeSent;
131 std::vector<int> indexToBeReceived;
132 operator bool() const {return on;}
DomainDecomposition()133 DomainDecomposition():
134 on(false), async(false)
135 {}
136 void enable(Communicator& c);
137 };
138
139 DomainDecomposition dd;
140 long int ddStep; //last step in which dd happened
141
142 void share(const std::set<AtomNumber>&);
143
144 public:
145
146 explicit Atoms(PlumedMain&plumed);
147 ~Atoms();
148
149 void init();
150
151 void share();
152 void shareAll();
153 void wait();
154 void updateForces();
155
156 void setRealPrecision(int);
157 int getRealPrecision()const;
158
159 void setTimeStep(void*);
160 double getTimeStep()const;
161
162 void setKbT(void*);
163 double getKbT()const;
164
165 void setNatoms(int);
166 int getNatoms()const;
167 int getNVirtualAtoms()const;
168
169 const long int& getDdStep()const;
170 const std::vector<int>& getGatindex()const;
171 const Pbc& getPbc()const;
172 void getLocalMasses(std::vector<double>&);
173 void getLocalPositions(std::vector<Vector>&);
174 void getLocalForces(std::vector<Vector>&);
175 void getLocalMDForces(std::vector<Vector>&);
176 const Tensor& getVirial()const;
177
setCollectEnergy(bool b)178 void setCollectEnergy(bool b) { collectEnergy=b; }
179
180 void setDomainDecomposition(Communicator&);
181 void setAtomsGatindex(int*,bool);
182 void setAtomsContiguous(int);
183 void setAtomsNlocal(int);
184
185 void startStep();
186 void setEnergy(void*);
187 void setBox(void*);
188 void setVirial(void*);
189 void setPositions(void*);
190 void setPositions(void*,int);
191 void setForces(void*);
192 void setForces(void*,int);
193 void setMasses(void*);
194 void setCharges(void*);
195 bool chargesWereSet() const ;
196 bool boxWasSet() const ;
197
198 void MD2double(const void*m,double&d)const;
199 void double2MD(const double&d,void*m)const;
200
201 void createFullList(int*);
202 void getFullList(int**);
203 void clearFullList();
204
205 void add(ActionAtomistic*);
206 void remove(ActionAtomistic*);
207
getEnergy()208 double getEnergy()const {plumed_assert(collectEnergy && energyHasBeenSet); return energy;}
209
isEnergyNeeded()210 bool isEnergyNeeded()const {return collectEnergy;}
211
setMDEnergyUnits(double d)212 void setMDEnergyUnits(double d) {MDUnits.setEnergy(d);}
setMDLengthUnits(double d)213 void setMDLengthUnits(double d) {MDUnits.setLength(d);}
setMDTimeUnits(double d)214 void setMDTimeUnits(double d) {MDUnits.setTime(d);}
setMDChargeUnits(double d)215 void setMDChargeUnits(double d) {MDUnits.setCharge(d);}
setMDMassUnits(double d)216 void setMDMassUnits(double d) {MDUnits.setMass(d);}
getMDUnits()217 const Units& getMDUnits() {return MDUnits;}
setUnits(const Units & u)218 void setUnits(const Units&u) {units=u;}
getUnits()219 const Units& getUnits() {return units;}
220 void updateUnits();
221
222 AtomNumber addVirtualAtom(ActionWithVirtualAtom*);
223 void removeVirtualAtom(ActionWithVirtualAtom*);
224 ActionWithVirtualAtom* getVirtualAtomsAction(AtomNumber)const;
225 bool isVirtualAtom(AtomNumber)const;
226 void insertGroup(const std::string&name,const std::vector<AtomNumber>&a);
227 void removeGroup(const std::string&name);
228 void writeBinary(std::ostream&)const;
229 void readBinary(std::istream&);
230 double getKBoltzmann()const;
231 double getMDKBoltzmann()const;
232 bool usingNaturalUnits()const;
setNaturalUnits(bool n)233 void setNaturalUnits(bool n) {naturalUnits=n;}
setMDNaturalUnits(bool n)234 void setMDNaturalUnits(bool n) {MDnaturalUnits=n;}
235
236 void setExtraCV(const std::string &name,void*p);
237 void setExtraCVForce(const std::string &name,void*p);
238 double getExtraCV(const std::string &name);
239 void updateExtraCVForce(const std::string &name,double f);
240 };
241
242 inline
getNatoms()243 int Atoms::getNatoms()const {
244 return natoms;
245 }
246
247 inline
getNVirtualAtoms()248 int Atoms::getNVirtualAtoms()const {
249 return virtualAtomsActions.size();
250 }
251
252 inline
getDdStep()253 const long int& Atoms::getDdStep()const {
254 return ddStep;
255 }
256
257 inline
getGatindex()258 const std::vector<int>& Atoms::getGatindex()const {
259 return gatindex;
260 }
261
262 inline
getPbc()263 const Pbc& Atoms::getPbc()const {
264 return pbc;
265 }
266
267 inline
isVirtualAtom(AtomNumber i)268 bool Atoms::isVirtualAtom(AtomNumber i)const {
269 return i.index()>=(unsigned) getNatoms();
270 }
271
272 inline
getVirtualAtomsAction(AtomNumber i)273 ActionWithVirtualAtom* Atoms::getVirtualAtomsAction(AtomNumber i)const {
274 return virtualAtomsActions[i.index()-getNatoms()];
275 }
276
277 inline
usingNaturalUnits()278 bool Atoms::usingNaturalUnits() const {
279 return naturalUnits || MDnaturalUnits;
280 }
281
282 inline
chargesWereSet()283 bool Atoms::chargesWereSet() const {
284 return chargesHaveBeenSet;
285 }
286
287 inline
boxWasSet()288 bool Atoms::boxWasSet() const {
289 return boxHasBeenSet;
290 }
291
292 inline
getVirial()293 const Tensor& Atoms::getVirial()const {
294 return virial;
295 }
296
297
298 }
299 #endif
300