Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2013-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_multicolvar_MultiColvarBase_h
23 : #define __PLUMED_multicolvar_MultiColvarBase_h
24 :
25 : #include "core/ActionAtomistic.h"
26 : #include "core/ActionWithValue.h"
27 : #include "tools/DynamicList.h"
28 : #include "tools/LinkCells.h"
29 : #include "vesselbase/StoreDataVessel.h"
30 : #include "vesselbase/ActionWithVessel.h"
31 : #include "CatomPack.h"
32 : #include <vector>
33 :
34 : namespace PLMD {
35 : namespace multicolvar {
36 :
37 : class AtomValuePack;
38 : class BridgedMultiColvarFunction;
39 : class ActionVolume;
40 :
41 : class MultiColvarBase :
42 : public ActionAtomistic,
43 : public ActionWithValue,
44 : public vesselbase::ActionWithVessel
45 : {
46 : friend class BridgedMultiColvarFunction;
47 : friend class VolumeGradientBase;
48 : friend class MultiColvarFilter;
49 : friend class AtomValuePack;
50 : private:
51 : /// Use periodic boundary conditions
52 : bool usepbc;
53 : /// The forces we are going to apply to things
54 : std::vector<double> forcesToApply;
55 : /// We use this to say that all the atoms in the third block should are in the tasks
56 : bool allthirdblockintasks;
57 : /// In certain cases we can make three atom link cells faster
58 : bool uselinkforthree;
59 : /// Number of atoms that are active on this step
60 : unsigned nactive_atoms;
61 : /// Stuff for link cells - this is used to make coordination number like variables faster
62 : LinkCells linkcells;
63 : /// Link cells for third block of atoms
64 : LinkCells threecells;
65 : /// Number of atoms that are being used for central atom position
66 : unsigned ncentral;
67 : /// Bool vector telling us which atoms are required to calculate central atom position
68 : std::vector<bool> use_for_central_atom;
69 : /// Vector of tempory holders for central atom values
70 : std::vector<CatomPack> my_tmp_capacks;
71 : /// 1/number of atoms involved in central atoms
72 : double numberForCentralAtom;
73 : /// Ensures that setup is only performed once per loop
74 : bool setup_completed;
75 : /// Ensures that retrieving of atoms is only done once per calculation loop
76 : bool atomsWereRetrieved;
77 : /// Add derivatives of center of mass position
78 : void addComDerivatives( const int& ival, const unsigned& iatom, const Vector& der, multicolvar::AtomValuePack& myatoms ) const ;
79 : protected:
80 : /// This is used to keep track of what is calculated where
81 : std::vector< std::pair<unsigned,unsigned> > atom_lab;
82 : /// The vessels in these multicolvars in which the data is stored
83 : std::vector<vesselbase::StoreDataVessel*> mybasedata;
84 : /// The multicolvars from which we construct these quantities
85 : std::vector<MultiColvarBase*> mybasemulticolvars;
86 : /// This remembers where the boundaries are for the tasks. It makes link cells work fast
87 : Matrix<std::pair<unsigned,unsigned> > bookeeping;
88 : /// Function that recursively checks if filters have been used in the input to a multicolvar
89 : /// we need this to ensure that setupLinkCells is run in calculate with some actions
90 : bool filtersUsedAsInput();
91 : /// This resizes the arrays that are used for link cell update
92 : void resizeBookeepingArray( const unsigned& num1, const unsigned& num2 );
93 : /// Are we doing sums of matrix rows
94 : bool matsums;
95 : /// Using the species keyword to read in atoms
96 : bool usespecies;
97 : /// Number of atoms in each block
98 : unsigned nblock;
99 : /// This is used when turning cvcodes into atom numbers
100 : std::vector<unsigned> decoder;
101 : /// Blocks of atom numbers
102 : std::vector< std::vector<unsigned> > ablocks;
103 : /// Add a task to the list of tasks
104 : void addTaskToList( const unsigned& taskCode );
105 : /// Finish setting up the multicolvar base
106 : void setupMultiColvarBase( const std::vector<AtomNumber>& atoms );
107 : /// Add some derivatives to a particular component of a particular atom
108 : void addAtomDerivatives( const int&, const unsigned&, const Vector&, multicolvar::AtomValuePack& ) const ;
109 : /// Add derivative of the input value
110 : void mergeInputDerivatives( const unsigned& ival, const unsigned& start, const unsigned& end, const unsigned& jatom,
111 : const std::vector<double>& der, MultiValue& myder, AtomValuePack& myatoms ) const ;
112 : /// This routine take the ith set of input derivatives and adds it to each of the (end-start) output derivatives
113 : /// In other words one set of derivatives comes in and end-start sets of derivatives come out
114 : void splitInputDerivatives( const unsigned& ival, const unsigned& start, const unsigned& end,
115 : const unsigned& jatom, const std::vector<double>& der,
116 : MultiValue& myder, AtomValuePack& myatoms ) const ;
117 : /// This is used to accumulate functions of the coordination sphere. Ensures weights are taken into account
118 : void accumulateSymmetryFunction( const int& ival, const unsigned& iatom, const double& val, const Vector& der, const Tensor& vir, multicolvar::AtomValuePack& myatoms ) const ;
119 : /// Set which atoms are to be used to calculate the central atom position
120 : void setAtomsForCentralAtom( const std::vector<bool>& catom_ind );
121 : /// Set the value of the cutoff for the link cells
122 : void setLinkCellCutoff( const double& lcut, double tcut=-1.0 );
123 : /// Setup the link cells and neighbour list stuff
124 : void setupActiveTaskSet( std::vector<unsigned>& active_tasks, const std::string& input_label );
125 : /// Setup link cells in order to make this calculation faster
126 : void setupLinkCells();
127 : /// Get the cutoff for the link cells
128 : double getLinkCellCutoff() const ;
129 : /// This does setup of link cell stuff that is specific to the non-use of the usespecies keyword
130 : void setupNonUseSpeciesLinkCells( const unsigned& );
131 : /// Get the separation between a pair of vectors
132 : Vector getSeparation( const Vector& vec1, const Vector& vec2 ) const ;
133 : /// This sets up the list of atoms that are involved in this colvar
134 : bool setupCurrentAtomList( const unsigned& taskCode, AtomValuePack& myatoms ) const ;
135 : /// Decode indices if there are 2 or 3 atoms involved
136 : void decodeIndexToAtoms( const unsigned& taskCode, std::vector<unsigned>& atoms ) const ;
137 : /// Read in some atoms
138 : bool parseMultiColvarAtomList(const std::string& key, const int& num, std::vector<AtomNumber>& t);
139 : /// Read in ATOMS keyword
140 : void readAtomsLikeKeyword( const std::string & key, const int& natoms, std::vector<AtomNumber>& all_atoms );
141 : /// Read in two groups of atoms and setup multicolvars to calculate
142 : void readTwoGroups( const std::string& key0, const std::string& key1, const std::string& key2, std::vector<AtomNumber>& all_atoms );
143 : /// Read in three groups of atoms
144 : void readGroupKeywords( const std::string& key0, const std::string& key1, const std::string& key2, const std::string& key3,
145 : const bool& no_third_dim_accum, const bool& symmetric, std::vector<AtomNumber>& all_atoms );
146 : /// Read in three groups of atoms and construct CVs involving at least one
147 : void readThreeGroups( const std::string& key1, const std::string& key2, const std::string& key3,
148 : const bool& allow2, const bool& no_third_dim_accum, std::vector<AtomNumber>& all_atoms );
149 : /// Build sets by taking one multicolvar from each base
150 : void buildSets();
151 : public:
152 : explicit MultiColvarBase(const ActionOptions&);
153 1312 : ~MultiColvarBase() {}
154 : static void registerKeywords( Keywords& keys );
155 : /// Turn on the derivatives
156 : virtual void turnOnDerivatives();
157 : /// Do we use pbc to calculate this quantity
158 : bool usesPbc() const ;
159 : /// Apply PBCs over a set of distance vectors
160 : void applyPbc(std::vector<Vector>& dlist, unsigned max_index=0) const;
161 : /// Is it safe to use multithreading
162 8874 : bool threadSafe() const { return !(mybasemulticolvars.size()>0); }
163 : /// Do some setup before the calculation
164 : void prepare();
165 : /// This is overwritten here in order to make sure that we do not retrieve atoms multiple times
166 : void retrieveAtoms();
167 : /// Do the calculation
168 : virtual void calculate();
169 : /// Calculate numerical derivatives
170 : virtual void calculateNumericalDerivatives( ActionWithValue* a=NULL );
171 : /// Perform one of the tasks
172 : virtual void performTask( const unsigned&, const unsigned&, MultiValue& ) const ;
173 : /// Update the active atoms
174 : virtual void updateActiveAtoms( AtomValuePack& myatoms ) const ;
175 : /// This gets the position of an atom for the link cell setup
176 : virtual Vector getPositionOfAtomForLinkCells( const unsigned& iatom ) const ;
177 : /// Get the absolute index of the central atom
178 : virtual AtomNumber getAbsoluteIndexOfCentralAtom( const unsigned& i ) const ;
179 : /// This is replaced once we have a function to calculate the cv
180 : virtual double compute( const unsigned& tindex, AtomValuePack& myatoms ) const=0;
181 : /// Apply the forces from this action
182 : virtual void apply();
183 : /// Get the number of derivatives for this action
184 : virtual unsigned getNumberOfDerivatives(); // N.B. This is replacing the virtual function in ActionWithValue
185 : /// Checks if an task is being performed at the present time
186 : virtual bool isCurrentlyActive( const unsigned& code );
187 : ///
188 : virtual void getCentralAtomPack( const unsigned& basn, const unsigned& curr, CatomPack& mypack);
189 : /// Get the index where the central atom is stored
190 : virtual Vector getCentralAtomPos( const unsigned& curr );
191 : /// You can use this to screen contributions that are very small so we can avoid expensive (and pointless) calculations
192 : virtual double calculateWeight( const unsigned& taskCode, const double& weight, AtomValuePack& myatoms ) const ;
193 : /// Is this a density?
194 1828438 : virtual bool isDensity() const { return false; }
195 : /// Is the iatom'th stored value currently active
196 : bool storedValueIsActive( const unsigned& iatom );
197 : /// This is true if multicolvar is calculating a vector or if the multicolvar is the density
198 0 : virtual bool hasDifferentiableOrientation() const { return false; }
199 : /// This makes sure we are not calculating the director when we do LocalAverage
200 0 : virtual void doNotCalculateDirector() {}
201 : /// Ensure that derivatives are only calculated when needed
202 : bool doNotCalculateDerivatives() const ;
203 : /// Get the icolv th base multicolvar
204 : MultiColvarBase* getBaseMultiColvar( const unsigned& icolv ) const ;
205 : /// Get the number of base multicolvars
206 : unsigned getNumberOfBaseMultiColvars() const ;
207 : /// Get an input data
208 : virtual void getInputData( const unsigned& ind, const bool& normed, const multicolvar::AtomValuePack& myatoms, std::vector<double>& orient ) const ;
209 : /// Retrieve the input derivatives
210 : virtual MultiValue& getInputDerivatives( const unsigned& iatom, const bool& normed, const multicolvar::AtomValuePack& myatoms ) const ;
211 : };
212 :
213 : inline
214 914028 : bool MultiColvarBase::isCurrentlyActive( const unsigned& code ) {
215 1085337 : if( setup_completed && atom_lab[code].first>0 ) {
216 17638 : unsigned mmc=atom_lab[code].first - 1;
217 35276 : return mybasedata[mmc]->storedValueIsActive( atom_lab[code].second );
218 : }
219 : return true;
220 : }
221 :
222 : inline
223 156 : AtomNumber MultiColvarBase::getAbsoluteIndexOfCentralAtom( const unsigned& iatom ) const {
224 : plumed_dbg_assert( iatom<atom_lab.size() );
225 312 : if( atom_lab[iatom].first>0 ) {
226 0 : unsigned mmc=atom_lab[iatom].first - 1;
227 0 : return mybasemulticolvars[mmc]->getAbsoluteIndexOfCentralAtom( atom_lab[iatom].second );
228 : }
229 : plumed_dbg_assert( usespecies );
230 312 : return ActionAtomistic::getAbsoluteIndex( atom_lab[getTaskCode(iatom)].second );
231 : }
232 :
233 : inline
234 405207244 : Vector MultiColvarBase::getPositionOfAtomForLinkCells( const unsigned& iatom ) const {
235 : plumed_dbg_assert( iatom<atom_lab.size() );
236 810414488 : if( atom_lab[iatom].first>0 ) {
237 2487901 : unsigned mmc=atom_lab[iatom].first - 1;
238 4975802 : return mybasemulticolvars[mmc]->getCentralAtomPos( atom_lab[iatom].second );
239 : }
240 805438686 : return ActionAtomistic::getPosition( atom_lab[iatom].second );
241 : }
242 :
243 : inline
244 8908998 : unsigned MultiColvarBase::getNumberOfDerivatives() {
245 8908998 : return 3*getNumberOfAtoms()+9;
246 : }
247 :
248 : inline
249 : bool MultiColvarBase::usesPbc() const {
250 234812 : return usepbc;
251 : }
252 :
253 : inline
254 65674493 : bool MultiColvarBase::doNotCalculateDerivatives() const {
255 65674493 : if( !dertime ) return true;
256 64987873 : return ActionWithValue::doNotCalculateDerivatives();
257 : }
258 :
259 : inline
260 : unsigned MultiColvarBase::getNumberOfBaseMultiColvars() const {
261 140 : return mybasemulticolvars.size();
262 : }
263 :
264 : inline
265 : MultiColvarBase* MultiColvarBase::getBaseMultiColvar( const unsigned& icolv ) const {
266 : plumed_dbg_assert( icolv<mybasemulticolvars.size() );
267 25007 : return mybasemulticolvars[icolv];
268 : }
269 :
270 : }
271 : }
272 :
273 : #endif
|
