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