1 /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 Copyright (c) 2012-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_vesselbase_ActionWithVessel_h
23 #define __PLUMED_vesselbase_ActionWithVessel_h
24
25 #include "core/ActionWithValue.h"
26 #include "core/ActionAtomistic.h"
27 #include "tools/Exception.h"
28 #include "tools/DynamicList.h"
29 #include "tools/MultiValue.h"
30 #include <vector>
31 #include "tools/ForwardDecl.h"
32
33 namespace PLMD {
34 class Value;
35 class Stopwatch;
36
37 namespace vesselbase {
38
39 class Vessel;
40 class BridgeVessel;
41 class StoreDataVessel;
42
43 /**
44 \ingroup MULTIINHERIT
45 This is used to create PLMD::Action objects that are computed by calculating the same function multiple
46 times. This is used in PLMD::MultiColvar.
47 */
48
49 class ActionWithVessel : public virtual Action {
50 friend class Vessel;
51 friend class ShortcutVessel;
52 friend class FunctionVessel;
53 friend class StoreDataVessel;
54 friend class BridgeVessel;
55 friend class ActionWithInputVessel;
56 friend class OrderingVessel;
57 private:
58 /// Do all calculations in serial
59 bool serial;
60 /// Lower memory requirements
61 bool lowmem;
62 /// Are we skipping the calculation of the derivatives
63 bool noderiv;
64 /// This tells plumed that this is used in a bridge
65 bool actionIsBridged;
66 /// The maximum number of derivatives we can use before we need to invoke lowmem
67 unsigned maxderivatives;
68 /// The tolerance on the accumulators
69 double tolerance;
70 /// Tolerance for quantities being put in neighbor lists
71 double nl_tolerance;
72 /// Pointers to the functions we are using on each value
73 std::vector<std::unique_ptr<Vessel>> functions;
74 /// Tempory storage for forces
75 std::vector<double> tmpforces;
76 /// Ths full list of tasks we have to perform
77 std::vector<unsigned> fullTaskList;
78 /// The current number of active tasks
79 unsigned nactive_tasks;
80 /// The indices of the tasks in the full list of tasks
81 std::vector<unsigned> indexOfTaskInFullList;
82 /// The list of currently active tasks
83 std::vector<unsigned> partialTaskList;
84 /// The list of atoms involved in derivatives (we keep a copy here to avoid resizing)
85 std::vector<unsigned> der_list;
86 /// The buffer that we use (we keep a copy here to avoid resizing)
87 std::vector<double> buffer;
88 /// Do we want to output information on the timings of different parts of the calculation
89 bool timers;
90 ForwardDecl<Stopwatch> stopwatch_fwd;
91 /// The stopwatch that times the different parts of the calculation
92 Stopwatch& stopwatch=*stopwatch_fwd;
93 /// These are used to minmise computational expense in complex functions
94 bool dertime_can_be_off;
95 protected:
96 /// This is also used to minimise computational expense in complex functions
97 bool dertime;
98 /// The terms in the series are locked
99 bool contributorsAreUnlocked;
100 /// Does the weight have derivatives
101 bool weightHasDerivatives;
102 /// This is used for numerical derivatives of bridge variables
103 unsigned bridgeVariable;
104 /// A pointer to the object that stores data
105 StoreDataVessel* mydata;
106 /// This list is used to update the neighbor list
107 std::vector<unsigned> taskFlags;
108 /// Add a vessel to the list of vessels
109 void addVessel( const std::string& name, const std::string& input, const int numlab=0 );
110 void addVessel( std::unique_ptr<Vessel> vv );
111 /// Add a bridging vessel to the list of vessels
112 BridgeVessel* addBridgingVessel( ActionWithVessel* tome );
113 /// Complete the setup of this object (this routine must be called after construction of ActionWithValue)
114 void readVesselKeywords();
115 /// Turn on the derivatives in the vessel
116 void needsDerivatives();
117 /// Return the value of the tolerance
118 double getTolerance() const ;
119 /// Return the value for the neighbor list tolerance
120 double getNLTolerance() const ;
121 /// Calculate the values of all the vessels
122 void runAllTasks();
123 /// Resize all the functions when the number of derivatives change
124 void resizeFunctions();
125 /// This loops over all the vessels calculating them and also
126 /// sets all the element derivatives equal to zero
127 void calculateAllVessels( const unsigned& taskCode, MultiValue& myvals, MultiValue& bvals, std::vector<double>& buffer, std::vector<unsigned>& der_list );
128 /// Retrieve the forces from all the vessels (used in apply)
129 bool getForcesFromVessels( std::vector<double>& forcesToApply );
130 /// Is the calculation being done in serial
131 bool serialCalculation() const;
132 /// Are we using low memory
133 bool usingLowMem() const ;
134 /// Set that we are using low memory
135 void setLowMemOption(const bool& );
136 /// Deactivate all the tasks in the task list
137 void deactivateAllTasks();
138 /// Get the size of the buffer
139 unsigned getSizeOfBuffer( unsigned& bufsize );
140 /// Add a task to the full list
141 void addTaskToList( const unsigned& taskCode );
142 public:
143 static void registerKeywords(Keywords& keys);
144 explicit ActionWithVessel(const ActionOptions&ao);
145 ~ActionWithVessel();
146 void lockContributors();
147 /// Get the number of tasks that are currently active
148 unsigned getCurrentNumberOfActiveTasks() const ;
149 /// Check whether or not a particular task is currently active
150 bool taskIsCurrentlyActive( const unsigned& index ) const ;
151 /// Are derivatives required for this quantity
152 bool derivativesAreRequired() const ;
153 /// Is this action thread safe
threadSafe()154 virtual bool threadSafe() const { return true; }
155 /// Finish running all the calculations
156 virtual void finishComputations( const std::vector<double>& buffer );
157 /// Are the base quantities periodic
158 virtual bool isPeriodic()=0;
159 /// What are the domains of the base quantities
160 virtual void retrieveDomain( std::string& min, std::string& max);
161 /// Get the number of derivatives for final calculated quantity
162 virtual unsigned getNumberOfDerivatives()=0;
163 /// Get the number of quantities that are calculated during each task
164 virtual unsigned getNumberOfQuantities() const ;
165 /// Get the number of vessels
166 unsigned getNumberOfVessels() const;
167 /// Get a pointer to the ith vessel
168 Vessel* getPntrToVessel( const unsigned& i );
169 /// Do any jobs that are required before the task list is undertaken
170 virtual void doJobsRequiredBeforeTaskList();
171 /// Get the full size of the taskList dynamic list
172 unsigned getFullNumberOfTasks() const ;
173 /// Get the position of the ith active task in the full list
174 unsigned getPositionInFullTaskList( const unsigned& ii ) const ;
175 /// Get the code for the ii th task in the list
176 unsigned getTaskCode( const unsigned& ii ) const ;
177 /// Get the ith of the currently active tasks
178 unsigned getActiveTask( const unsigned& ii ) const ;
179 /// Calculate one of the functions in the distribution
180 virtual void performTask( const unsigned&, const unsigned&, MultiValue& ) const=0;
181 /// Do the task if we have a bridge
182 virtual void transformBridgedDerivatives( const unsigned& current, MultiValue& invals, MultiValue& outvals ) const;
183 /// Ensure that data required in other vessels is stored
184 StoreDataVessel* buildDataStashes( ActionWithVessel* actionThatUses );
185 /// Apply forces from bridge vessel - this is rarely used - currently only in ActionVolume
applyBridgeForces(const std::vector<double> & bb)186 virtual void applyBridgeForces( const std::vector<double>& bb ) { plumed_error(); }
187 /// These are overwritten in MultiColvarFunction
188 // virtual void activateIndexes( const unsigned&, const unsigned&, const std::vector<unsigned>& ){}
189 /// Return a particular named vessel
190 Vessel* getVesselWithName( const std::string& mynam );
191 /// Does the weight have derivatives
192 bool weightWithDerivatives() const ;
193 /// Return the position in the current task list
194 unsigned getPositionInCurrentTaskList( const unsigned& myind ) const ;
195 /// These normalizes vectors and is used in StoreDataVessel
normalizeVector(std::vector<double> & vals)196 virtual void normalizeVector( std::vector<double>& vals ) const { plumed_error(); }
normalizeVectorDerivatives(MultiValue & myvals)197 virtual void normalizeVectorDerivatives( MultiValue& myvals ) const { plumed_error(); }
198 };
199
200 inline
getTolerance()201 double ActionWithVessel::getTolerance() const {
202 return tolerance;
203 }
204
205 inline
getNLTolerance()206 double ActionWithVessel::getNLTolerance() const {
207 return nl_tolerance;
208 }
209
210 inline
getNumberOfVessels()211 unsigned ActionWithVessel::getNumberOfVessels() const {
212 return functions.size();
213 }
214
215 inline
getNumberOfQuantities()216 unsigned ActionWithVessel::getNumberOfQuantities() const {
217 return 2;
218 }
219
220 inline
getPntrToVessel(const unsigned & i)221 Vessel* ActionWithVessel::getPntrToVessel( const unsigned& i ) {
222 plumed_dbg_assert( i<functions.size() );
223 return functions[i].get();
224 }
225
226 inline
getFullNumberOfTasks()227 unsigned ActionWithVessel::getFullNumberOfTasks() const {
228 return fullTaskList.size();
229 }
230
231 inline
getTaskCode(const unsigned & ii)232 unsigned ActionWithVessel::getTaskCode( const unsigned& ii ) const {
233 plumed_dbg_assert( ii<fullTaskList.size() );
234 return fullTaskList[ii];
235 }
236
237 inline
getCurrentNumberOfActiveTasks()238 unsigned ActionWithVessel::getCurrentNumberOfActiveTasks() const {
239 return nactive_tasks;
240 }
241
242 inline
getActiveTask(const unsigned & ii)243 unsigned ActionWithVessel::getActiveTask( const unsigned& ii ) const {
244 plumed_dbg_assert( ii<nactive_tasks );
245 return partialTaskList[ii];
246 }
247
248 inline
getPositionInFullTaskList(const unsigned & ii)249 unsigned ActionWithVessel::getPositionInFullTaskList( const unsigned& ii ) const {
250 plumed_dbg_assert( ii<nactive_tasks );
251 return indexOfTaskInFullList[ii];
252 }
253
254 inline
serialCalculation()255 bool ActionWithVessel::serialCalculation() const {
256 return serial;
257 }
258
259 inline
usingLowMem()260 bool ActionWithVessel::usingLowMem() const {
261 return lowmem;
262 }
263
264 inline
setLowMemOption(const bool & l)265 void ActionWithVessel::setLowMemOption(const bool& l) {
266 lowmem=l;
267 }
268
269 inline
derivativesAreRequired()270 bool ActionWithVessel::derivativesAreRequired() const {
271 return !noderiv;
272 }
273
274 inline
weightWithDerivatives()275 bool ActionWithVessel::weightWithDerivatives() const {
276 return weightHasDerivatives;
277 }
278
279 inline
getPositionInCurrentTaskList(const unsigned & myind)280 unsigned ActionWithVessel::getPositionInCurrentTaskList( const unsigned& myind ) const {
281 if( nactive_tasks==fullTaskList.size() ) return myind;
282
283 for(unsigned i=0; i<nactive_tasks; ++i) {
284 if( myind==indexOfTaskInFullList[i] ) return i;
285 }
286 plumed_merror("requested task is not active");
287 }
288
289 }
290 }
291 #endif
292