Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-2019 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_ActionAtomistic_h
23 : #define __PLUMED_core_ActionAtomistic_h
24 :
25 : #include "Action.h"
26 : #include "tools/Tensor.h"
27 : #include "Atoms.h"
28 : #include "tools/Pbc.h"
29 : #include <vector>
30 : #include <set>
31 :
32 : namespace PLMD {
33 :
34 : class Pbc;
35 : class PDB;
36 :
37 : /// \ingroup MULTIINHERIT
38 : /// Action used to create objects that access the positions of the atoms from the MD code
39 : class ActionAtomistic :
40 : virtual public Action
41 : {
42 :
43 : std::vector<AtomNumber> indexes; // the set of needed atoms
44 : /// unique should be an ordered set since we later create a vector containing the corresponding indexes
45 : std::set<AtomNumber> unique;
46 : /// unique_local should be an ordered set since we later create a vector containing the corresponding indexes
47 : std::set<AtomNumber> unique_local;
48 : std::vector<Vector> positions; // positions of the needed atoms
49 : double energy;
50 : Pbc& pbc;
51 : Tensor virial;
52 : std::vector<double> masses;
53 : bool chargesWereSet;
54 : std::vector<double> charges;
55 :
56 : std::vector<Vector> forces; // forces on the needed atoms
57 : double forceOnEnergy;
58 :
59 : bool lockRequestAtoms; // forbid changes to request atoms
60 :
61 : bool donotretrieve;
62 : bool donotforce;
63 :
64 : protected:
65 : Atoms& atoms;
66 :
67 : public:
68 : /// Request an array of atoms.
69 : /// This method is used to ask for a list of atoms. Atoms
70 : /// should be asked for by number. If this routine is called
71 : /// during the simulation, atoms will be available at the next step
72 : /// MAYBE WE HAVE TO FIND SOMETHING MORE CLEAR FOR DYNAMIC
73 : /// LISTS OF ATOMS
74 : void requestAtoms(const std::vector<AtomNumber> & a);
75 : /// Get position of i-th atom (access by relative index)
76 : const Vector & getPosition(int)const;
77 : /// Get position of i-th atom (access by absolute AtomNumber).
78 : /// With direct access to the global atom array
79 : const Vector & getPosition(AtomNumber)const;
80 : /// Get modifiable position of i-th atom (access by absolute AtomNumber).
81 : /// Should be used by action that need to modify the stored atomic coordinates
82 : Vector & modifyPosition(AtomNumber);
83 : /// Get total number of atoms, including virtual ones.
84 : /// Can be used to make a loop on modifyPosition or getPosition(AtomNumber)
85 : unsigned getTotAtoms()const;
86 : /// Get modifiable force of i-th atom (access by absolute AtomNumber).
87 : /// \warning Should be used by action that need to modify the stored atomic forces.
88 : /// This should be used with great care since the plumed core does
89 : /// not usually keep all these forces up to date. In particular,
90 : /// if an action require this, one should during constructor
91 : /// call allowToAccessGlobalForces().
92 : /// Notice that for efficiency reason plumed does not check if this is done!
93 : Vector & modifyGlobalForce(AtomNumber);
94 : /// Get modifiable virial
95 : /// Should be used by action that need to modify the stored virial
96 : Tensor & modifyGlobalVirial();
97 : /// Get modifiable PBC
98 : /// Should be used by action that need to modify the stored box
99 : Pbc & modifyGlobalPbc();
100 : /// Get box shape
101 : const Tensor & getBox()const;
102 : /// Get the array of all positions
103 : const std::vector<Vector> & getPositions()const;
104 : /// Get energy
105 : const double & getEnergy()const;
106 : /// Get mass of i-th atom
107 : double getMass(int i)const;
108 : /// Get charge of i-th atom
109 : double getCharge(int i)const;
110 : /// Get a reference to forces array
111 : std::vector<Vector> & modifyForces();
112 : /// Get a reference to virial array
113 : Tensor & modifyVirial();
114 : /// Get a reference to force on energy
115 : double & modifyForceOnEnergy();
116 : /// Get number of available atoms
117 4185379809 : unsigned getNumberOfAtoms()const {return indexes.size();}
118 : /// Compute the pbc distance between two positions
119 : Vector pbcDistance(const Vector&,const Vector&)const;
120 : /// Applies PBCs to a seriens of positions or distances
121 : void pbcApply(std::vector<Vector>& dlist, unsigned max_index=0) const;
122 : /// Get the vector of absolute indexes
123 : virtual const std::vector<AtomNumber> & getAbsoluteIndexes()const;
124 : /// Get the absolute index of an atom
125 : AtomNumber getAbsoluteIndex(int i)const;
126 : /// Parse a list of atoms without a numbered keyword
127 : void parseAtomList(const std::string&key,std::vector<AtomNumber> &t);
128 : /// Parse an list of atom with a numbred keyword
129 : void parseAtomList(const std::string&key,const int num, std::vector<AtomNumber> &t);
130 : /// Convert a set of read in strings into an atom list (this is used in parseAtomList)
131 : void interpretAtomList( std::vector<std::string>& strings, std::vector<AtomNumber> &t);
132 : /// Change the box shape
133 : void changeBox( const Tensor& newbox );
134 : /// Get reference to Pbc
135 : const Pbc & getPbc() const;
136 : /// Add the forces to the atoms
137 : void setForcesOnAtoms( const std::vector<double>& forcesToApply, unsigned ind=0 );
138 : /// Skip atom retrieval - use with care.
139 : /// If this function is called during initialization, then atoms are
140 : /// not going to be retrieved. Can be used for optimization. Notice that
141 : /// calling getPosition(int) in an Action where DoNotRetrieve() was called might
142 : /// lead to undefined behavior.
143 25 : void doNotRetrieve() {donotretrieve=true;}
144 : /// Skip atom forces - use with care.
145 : /// If this function is called during initialization, then forces are
146 : /// not going to be propagated. Can be used for optimization.
147 17 : void doNotForce() {donotforce=true;}
148 : /// Make atoms whole, assuming they are in the proper order
149 : void makeWhole();
150 : /// Allow calls to modifyGlobalForce()
151 8 : void allowToAccessGlobalForces() {atoms.zeroallforces=true;}
152 : /// updates local unique atoms
153 : void updateUniqueLocal();
154 : public:
155 :
156 : // virtual functions:
157 :
158 : explicit ActionAtomistic(const ActionOptions&ao);
159 : ~ActionAtomistic();
160 :
161 : static void registerKeywords( Keywords& keys );
162 :
163 : void clearOutputForces();
164 :
165 : /// N.B. only pass an ActionWithValue to this routine if you know exactly what you
166 : /// are doing. The default will be correct for the vast majority of cases
167 : virtual void calculateNumericalDerivatives( ActionWithValue* a=NULL );
168 : /// Numerical derivative routine to use when using Actions that inherit from BOTH
169 : /// ActionWithArguments and ActionAtomistic
170 : void calculateAtomicNumericalDerivatives( ActionWithValue* a, const unsigned& startnum );
171 :
172 : virtual void retrieveAtoms();
173 : void applyForces();
174 : void lockRequests();
175 : void unlockRequests();
176 : const std::set<AtomNumber> & getUnique()const;
177 : const std::set<AtomNumber> & getUniqueLocal()const;
178 : /// Read in an input file containing atom positions and calculate the action for the atomic
179 : /// configuration therin
180 : void readAtomsFromPDB( const PDB& pdb );
181 : };
182 :
183 : inline
184 : const Vector & ActionAtomistic::getPosition(int i)const {
185 972797160 : return positions[i];
186 : }
187 :
188 : inline
189 : const Vector & ActionAtomistic::getPosition(AtomNumber i)const {
190 11478 : return atoms.positions[i.index()];
191 : }
192 :
193 : inline
194 : Vector & ActionAtomistic::modifyPosition(AtomNumber i) {
195 29186 : return atoms.positions[i.index()];
196 : }
197 :
198 : inline
199 : Vector & ActionAtomistic::modifyGlobalForce(AtomNumber i) {
200 14638 : return atoms.forces[i.index()];
201 : }
202 :
203 : inline
204 : Tensor & ActionAtomistic::modifyGlobalVirial() {
205 113 : return atoms.virial;
206 : }
207 :
208 : inline
209 : double ActionAtomistic::getMass(int i)const {
210 91026 : return masses[i];
211 : }
212 :
213 : inline
214 19683 : double ActionAtomistic::getCharge(int i) const {
215 19683 : if( !chargesWereSet ) error("charges were not passed to plumed");
216 39366 : return charges[i];
217 : }
218 :
219 : inline
220 1234 : const std::vector<AtomNumber> & ActionAtomistic::getAbsoluteIndexes()const {
221 1337 : return indexes;
222 : }
223 :
224 : inline
225 : AtomNumber ActionAtomistic::getAbsoluteIndex(int i)const {
226 260418353 : return indexes[i];
227 : }
228 :
229 : inline
230 : const std::vector<Vector> & ActionAtomistic::getPositions()const {
231 550286 : return positions;
232 : }
233 :
234 : inline
235 : const double & ActionAtomistic::getEnergy()const {
236 3787 : return energy;
237 : }
238 :
239 : inline
240 : const Tensor & ActionAtomistic::getBox()const {
241 24523 : return pbc.getBox();
242 : }
243 :
244 : inline
245 : std::vector<Vector> & ActionAtomistic::modifyForces() {
246 77295 : return forces;
247 : }
248 :
249 : inline
250 : Tensor & ActionAtomistic::modifyVirial() {
251 80640 : return virial;
252 : }
253 :
254 : inline
255 : double & ActionAtomistic::modifyForceOnEnergy() {
256 : return forceOnEnergy;
257 : }
258 :
259 : inline
260 : const Pbc & ActionAtomistic::getPbc() const {
261 310284 : return pbc;
262 : }
263 :
264 : inline
265 86038 : void ActionAtomistic::lockRequests() {
266 93190 : lockRequestAtoms=true;
267 86038 : }
268 :
269 : inline
270 86038 : void ActionAtomistic::unlockRequests() {
271 93190 : lockRequestAtoms=false;
272 86038 : }
273 :
274 : inline
275 : const std::set<AtomNumber> & ActionAtomistic::getUnique()const {
276 : return unique;
277 : }
278 :
279 : inline
280 : const std::set<AtomNumber> & ActionAtomistic::getUniqueLocal()const {
281 : return unique_local;
282 : }
283 :
284 : inline
285 : unsigned ActionAtomistic::getTotAtoms()const {
286 57662 : return atoms.positions.size();
287 : }
288 :
289 : inline
290 : Pbc & ActionAtomistic::modifyGlobalPbc() {
291 65 : return atoms.pbc;
292 : }
293 :
294 :
295 :
296 : }
297 :
298 : #endif
|