Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2012-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 : #include "core/ActionShortcut.h"
23 : #include "core/ActionRegister.h"
24 : #include "MultiColvarShortcuts.h"
25 : #include <string>
26 : #include <cmath>
27 :
28 : //+PLUMEDOC MCOLVAR DISTANCES
29 : /*
30 : Calculate the distances between multiple piars of atoms
31 :
32 : __This shortcut action allows you to calculate function of the distribution of distances and reproduces the syntax in older PLUMED versions.
33 : We think that the syntax in newer versions is more flexible and easier to use and that you should try to avoid using this old syntax. If you look at the example inputs below you can
34 : see how the new syntax can be used to replace this old one. We would strongly encourage you to use the newer syntax.__
35 :
36 : The following input tells plumed to calculate the distances between atoms 3 and 5 and
37 : between atoms 1 and 2 and to print the shorter of these two distances.
38 :
39 : ```plumed
40 : d1: DISTANCES ATOMS1=3,5 ATOMS2=1,2 LOWEST
41 : PRINT ARG=d1.lowest
42 : ```
43 :
44 : The following input tells plumed to calculate the distances between atoms 3 and 5 and between atoms 1 and 2
45 : and then to calculate the number of these distances that are less than 0.1 nm. The number of distances
46 : less than 0.1nm is then printed to a file.
47 :
48 : ```plumed
49 : d1: DISTANCES ATOMS1=3,5 ATOMS2=1,2 LESS_THAN={RATIONAL R_0=0.1}
50 : PRINT ARG=d1.lessthan
51 : ```
52 :
53 : The following input tells plumed to calculate all 3 distances between atoms 1, 2 and 3 (i.e. the distances between atoms
54 : 1 and 2, atoms 1 and 3 and atoms 2 and 3). The average of these distances is then calculated in two ways.
55 :
56 : ```plumed
57 : d1: DISTANCES GROUP=1-3 MEAN SUM
58 : d1mean: CUSTOM ARG=d1.sum FUNC=x/3 PERIODIC=NO
59 : PRINT ARG=d1.mean,d1mean
60 : ```
61 :
62 : The following input tells plumed to calculate all the distances between the atoms in GROUPA and the atoms in GROUPB.
63 : In other words the distances between atoms 1 and 2 and the distance between atoms 1 and 3. The number of distances
64 : more than 0.1 is then printed to a file.
65 :
66 : ```plumed
67 : d1: DISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
68 : PRINT ARG=d1.morethan
69 : ```
70 :
71 : ## Calculating minimum distances
72 :
73 : There are various way for calculating the minimum distance between two groups of atoms. The simplest one
74 : is the one we saw above; namely:
75 :
76 : ```plumed
77 : d1: DISTANCES GROUPA=1-10 GROUPB=11-20 LOWEST
78 : PRINT ARG=d1.lowest
79 : ```
80 :
81 : A disadvantage of this approach is that the derivatives of the resulting value are not continuous. To be clear,
82 : the value of the shortest distance is continuous when you use this approach so you will likely be fine if you use
83 : the quantity output by the command above as an input in a bias. However, if you would like to use a quantity with
84 : continuous derivatives for the minimum distance you can use the following input:
85 :
86 : ```plumed
87 : d1: DISTANCES GROUPA=1-10 GROUPB=11-20 MIN={BETA=500.}
88 : PRINT ARG=d1.min FILE=colvar STRIDE=10
89 : ```
90 :
91 : In order to ensure that the minimum value has continuous derivatives we use the following function:
92 :
93 : $$
94 : d_{min} = \frac{\beta}{ \log \sum_i \exp\left( \frac{\beta}{d_i} \right) }
95 : $$
96 :
97 : where $\beta$ is a user specified parameter.
98 :
99 : This input is used rather than a separate MINDIST colvar so that the same routine and the same input style can be
100 : used to calculate minimum coordination numbers (see [COORDINATIONNUMBER](COORDINATIONNUMBER.md)), minimum
101 : angles (see [ANGLES](ANGLES.md)) and many other variables. If you expand the input above you will notice that
102 : the DISTANCES command is a shortcut. By expanding this shortcut you can see how PLUMED calculates the minimum distance.
103 :
104 : An alternative implementation of the minimum distance is provided through the ALT_MIN keyword that is demonstrated below:
105 :
106 : ```plumed
107 : d1: DISTANCES GROUPA=1-10 GROUPB=11-20 ALT_MIN={BETA=500.}
108 : PRINT ARG=d1.altmin FILE=colvar STRIDE=10
109 : ```
110 :
111 : The minimum value output by this command also has continuous derivatives and is evaluated using the following function:
112 :
113 : $$
114 : d_{min} = -\frac{1}{\beta} \log \left( \sum_i \exp(-\beta d_i) \right)
115 : $$
116 :
117 : ## Calculating maximum distances
118 :
119 : There are two ways for calculating the maximum distance between two groups of atom. The simplest one is:
120 :
121 : ```plumed
122 : d1: DISTANCES GROUPA=1-10 GROUPB=11-20 HIGHEST
123 : PRINT ARG=d1.highest
124 : ```
125 :
126 : As with LOWEST above a disadvantage of this approach is that the derivatives of the resulting value are not continuous.
127 : If you would like to use a quantity with continuous derivatives for the maximum distance you can use the following input:
128 :
129 : ```plumed
130 : d1: DISTANCES GROUPA=1-10 GROUPB=11-20 MAX={BETA=500.}
131 : PRINT ARG=d1.max
132 : ```
133 :
134 : The maximum value output by this command is evaluated using the following function:
135 :
136 : $$
137 : d_{max} = \beta \log \left( \sum_i \exp\left( \frac{x}{\beta} \right) \right)
138 : $$
139 :
140 : ## Other parts of the old multicolvar syntax
141 :
142 : The ways of calculating the minimum and maximum distance described above are part of plumed 2's multicolvar functionality. In older
143 : versions of PLUMED these special actions
144 : allowed you to calculate multiple functions of a distribution of simple collective variables. As an example you
145 : can calculate the number of distances less than 1.0, the minimum distance, the number of distances more than
146 : 2.0, the number of distances between 1.0 and 2.0 and histogram that esimates the distribution of distances by using the following command:
147 :
148 : ```plumed
149 : d1: DISTANCES ...
150 : GROUPA=1-10 GROUPB=11-20
151 : LESS_THAN={RATIONAL R_0=1.0}
152 : MORE_THAN={RATIONAL R_0=2.0}
153 : BETWEEN={GAUSSIAN LOWER=1.0 UPPER=2.0}
154 : HISTOGRAM={GAUSSIAN LOWER=0.0 UPPER=4.0 NBINS=10}
155 : MIN={BETA=500.}
156 : ...
157 : PRINT ARG=d1.* FILE=colvar STRIDE=10
158 : ```
159 :
160 : A calculation performed this way is fast because the expensive part of the calculation - the calculation of all the distances - is only done once per step.
161 : If you expand the shortcut in the input above you can see how this result is achieved using a combination of many actions. We recommend you use these combinations of
162 : actions rather than the shortcut when using newer versions of PLUMED as you will have a more control over what you are doing. This newer input syntax is easier to
163 : understand especially when you start evaluating more complicated CVs.
164 :
165 :
166 : */
167 : //+ENDPLUMEDOC
168 :
169 : //+PLUMEDOC MCOLVAR XDISTANCES
170 : /*
171 : Calculate the x components of the vectors connecting one or many pairs of atoms.
172 :
173 : __This shortcut action allows you to calculate functions of the distribution of the x components of the vectors connecting pairs of atoms and
174 : reproduces the syntax in older PLUMED versions.
175 : If you look at the example inputs below you can
176 : see how the new syntax operates. We would strongly encourage you to use the newer syntax as it offers greater flexibility.__
177 :
178 : The following input tells plumed to calculate the x-component of the vector connecting atom 3 to atom 5 and
179 : the x-component of the vector connecting atom 1 to atom 2. The minimum of these two quantities is then
180 : printed
181 :
182 : ```plumed
183 : d1: XDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1}
184 : PRINT ARG=d1.min
185 : ```
186 :
187 : The following input tells plumed to calculate the x-component of the vector connecting atom 3 to atom 5 and
188 : the x-component of the vector connecting atom 1 to atom 2. The number of values that are
189 : less than 0.1nm is then printed to a file.
190 :
191 : ```plumed
192 : d1: XDISTANCES ATOMS1=3,5 ATOMS2=1,2 LESS_THAN={RATIONAL R_0=0.1}
193 : PRINT ARG=d1.lessthan
194 : ```
195 :
196 : The following input tells plumed to calculate the x-components of all the distinct vectors that can be created
197 : between atoms 1, 2 and 3 (i.e. the vectors between atoms 1 and 2, atoms 1 and 3 and atoms 2 and 3).
198 : The average of these quantities is then calculated.
199 :
200 : ```plumed
201 : d1: XDISTANCES GROUP=1-3 MEAN
202 : PRINT ARG=d1.mean
203 : ```
204 :
205 : The following input tells plumed to calculate all the vectors connecting the atoms in GROUPA to the atoms in GROUPB.
206 : In other words the vector between atoms 1 and 2 and the vector between atoms 1 and 3. The number of values
207 : more than 0.1 is then printed to a file.
208 :
209 : ```plumed
210 : d1: XDISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
211 : PRINT ARG=d1.morethan
212 : ```
213 :
214 : */
215 : //+ENDPLUMEDOC
216 :
217 :
218 : //+PLUMEDOC MCOLVAR YDISTANCES
219 : /*
220 : Calculate the y components of the vectors connecting one or many pairs of atoms.
221 :
222 : __This shortcut action allows you to calculate functions of the distribution of the y components of the vectors connecting pairs of atoms and
223 : reproduces the syntax in older PLUMED versions.
224 : If you look at the example inputs below you can
225 : see how the new syntax operates. We would strongly encourage you to use the newer syntax as it offers greater flexibility.__
226 :
227 : The following input tells plumed to calculate the y-component of the vector connecting atom 3 to atom 5 and
228 : the y-component of the vector connecting atom 1 to atom 2. The minimum of these two quantities is then
229 : printed
230 :
231 : ```plumed
232 : d1: YDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1}
233 : PRINT ARG=d1.min
234 : ```
235 :
236 : The following input tells plumed to calculate the y-component of the vector connecting atom 3 to atom 5 and
237 : the y-component of the vector connecting atom 1 to atom 2. The number of values that are
238 : less than 0.1nm is then printed to a file.
239 :
240 : ```plumed
241 : d1: YDISTANCES ATOMS1=3,5 ATOMS2=1,2 LESS_THAN={RATIONAL R_0=0.1}
242 : PRINT ARG=d1.lessthan
243 : ```
244 :
245 : The following input tells plumed to calculate the y-components of all the distinct vectors that can be created
246 : between atoms 1, 2 and 3 (i.e. the vectors between atoms 1 and 2, atoms 1 and 3 and atoms 2 and 3).
247 : The average of these quantities is then calculated.
248 :
249 : ```plumed
250 : d1: YDISTANCES GROUP=1-3 MEAN
251 : PRINT ARG=d1.mean
252 : ```
253 :
254 : The following input tells plumed to calculate all the vectors connecting the atoms in GROUPA to the atoms in GROUPB.
255 : In other words the vector between atoms 1 and 2 and the vector between atoms 1 and 3. The number of values
256 : more than 0.1 is then printed to a file.
257 :
258 : ```plumed
259 : d1: YDISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
260 : PRINT ARG=d1.morethan
261 : ```
262 :
263 : */
264 : //+ENDPLUMEDOC
265 :
266 : //+PLUMEDOC MCOLVAR ZDISTANCES
267 : /*
268 : Calculate the z components of the vectors connecting one or many pairs of atoms.
269 :
270 : __This shortcut action allows you to calculate functions of the distribution of the z components of the vectors connecting pairs of atoms and
271 : reproduces the syntax in older PLUMED versions.
272 : If you look at the example inputs below you can
273 : see how the new syntax operates. We would strongly encourage you to use the newer syntax as it offers greater flexibility.__
274 :
275 : The following input tells plumed to calculate the z-component of the vector connecting atom 3 to atom 5 and
276 : the z-component of the vector connecting atom 1 to atom 2. The minimum of these two quantities is then
277 : printed
278 :
279 : ```plumed
280 : d1: ZDISTANCES ATOMS1=3,5 ATOMS2=1,2 MIN={BETA=0.1}
281 : PRINT ARG=d1.min
282 : ```
283 :
284 : The following input tells plumed to calculate the z-component of the vector connecting atom 3 to atom 5 and
285 : the z-component of the vector connecting atom 1 to atom 2. The number of values that are
286 : less than 0.1nm is then printed to a file.
287 :
288 : ```plumed
289 : d1: ZDISTANCES ATOMS1=3,5 ATOMS2=1,2 LESS_THAN={RATIONAL R_0=0.1}
290 : PRINT ARG=d1.lessthan
291 : ```
292 :
293 : The following input tells plumed to calculate the z-components of all the distinct vectors that can be created
294 : between atoms 1, 2 and 3 (i.e. the vectors between atoms 1 and 2, atoms 1 and 3 and atoms 2 and 3).
295 : The average of these quantities is then calculated.
296 :
297 : ```plumed
298 : d1: ZDISTANCES GROUP=1-3 MEAN
299 : PRINT ARG=d1.mean
300 : ```
301 :
302 : The following input tells plumed to calculate all the vectors connecting the atoms in GROUPA to the atoms in GROUPB.
303 : In other words the vector between atoms 1 and 2 and the vector between atoms 1 and 3. The number of values
304 : more than 0.1 is then printed to a file.
305 :
306 : ```plumed
307 : d1: ZDISTANCES GROUPA=1 GROUPB=2,3 MORE_THAN={RATIONAL R_0=0.1}
308 : PRINT ARG=d1.morethan
309 : ```
310 :
311 : */
312 : //+ENDPLUMEDOC
313 :
314 :
315 : namespace PLMD {
316 : namespace multicolvar {
317 :
318 : class Distances : public ActionShortcut {
319 : public:
320 : static void registerKeywords(Keywords& keys);
321 : explicit Distances(const ActionOptions&);
322 : };
323 :
324 : PLUMED_REGISTER_ACTION(Distances,"DISTANCES")
325 : PLUMED_REGISTER_ACTION(Distances,"XDISTANCES")
326 : PLUMED_REGISTER_ACTION(Distances,"YDISTANCES")
327 : PLUMED_REGISTER_ACTION(Distances,"ZDISTANCES")
328 :
329 136 : void Distances::registerKeywords(Keywords& keys) {
330 136 : ActionShortcut::registerKeywords( keys );
331 272 : keys.setDeprecated("DISTANCE");
332 136 : keys.add("atoms-1","GROUP","Calculate the distance between each distinct pair of atoms in the group");
333 136 : keys.add("atoms-2","GROUPA","Calculate the distances between all the atoms in GROUPA and all "
334 : "the atoms in GROUPB. This must be used in conjunction with GROUPB.");
335 136 : keys.add("atoms-2","GROUPB","Calculate the distances between all the atoms in GROUPA and all the atoms "
336 : "in GROUPB. This must be used in conjunction with GROUPA.");
337 136 : keys.add("numbered","ATOMS","the pairs of atoms that you would like to calculate the angles for");
338 136 : keys.addFlag("NOPBC",false,"ignore the periodic boundary conditions when calculating distances");
339 136 : keys.addFlag("COMPONENTS",false,"calculate the x, y and z components of the distance separately and store them as label.x, label.y and label.z");
340 136 : keys.addFlag("SCALED_COMPONENTS",false,"calculate the a, b and c scaled components of the distance separately and store them as label.a, label.b and label.c");
341 136 : keys.addFlag("SQUARED",false,"calculate the x, y and z components of the distance separately and store them as label.x, label.y and label.z");
342 136 : keys.addFlag("USEGPU",false,"run this calculation on the GPU");
343 272 : keys.addDeprecatedFlag("LOWMEM","");
344 272 : keys.reset_style("ATOMS","atoms");
345 136 : MultiColvarShortcuts::shortcutKeywords( keys );
346 136 : keys.add("atoms","ORIGIN","calculate the distance of all the atoms specified using the ATOMS keyword from this point");
347 136 : keys.add("numbered","LOCATION","the location at which the CV is assumed to be in space");
348 272 : keys.reset_style("LOCATION","atoms");
349 272 : keys.setValueDescription("vector","the DISTANCES between the each pair of atoms that were specified");
350 272 : keys.addOutputComponent("x","COMPONENTS","vector","the x-components of the distance vectors");
351 272 : keys.addOutputComponent("y","COMPONENTS","vector","the y-components of the distance vectors");
352 272 : keys.addOutputComponent("z","COMPONENTS","vector","the z-components of the distance vectors");
353 136 : keys.needsAction("GROUP");
354 136 : keys.needsAction("DISTANCE");
355 136 : keys.needsAction("CENTER");
356 136 : }
357 :
358 71 : Distances::Distances(const ActionOptions& ao):
359 : Action(ao),
360 71 : ActionShortcut(ao) {
361 : // Create distances
362 : bool lowmem;
363 71 : parseFlag("LOWMEM",lowmem);
364 71 : if( lowmem ) {
365 0 : warning("LOWMEM flag is deprecated and is no longer required for this action");
366 : }
367 71 : std::string dline = getShortcutLabel() + ": DISTANCE";
368 : bool nopbc;
369 71 : parseFlag("NOPBC",nopbc);
370 71 : if( nopbc ) {
371 : dline += " NOPBC";
372 : }
373 71 : if( getName()=="DISTANCES" ) {
374 : bool comp;
375 68 : parseFlag("COMPONENTS",comp);
376 68 : if( comp ) {
377 : dline += " COMPONENTS";
378 : }
379 : bool scomp;
380 68 : parseFlag("SCALED_COMPONENTS",scomp);
381 68 : if( scomp ) {
382 : dline += " SCALED_COMPONENTS";
383 : }
384 : bool sq;
385 68 : parseFlag("SQUARED",sq);
386 68 : if( sq ) {
387 : dline += " SQUARED";
388 : }
389 : bool gpu;
390 73 : parseFlag("USEGPU",gpu);
391 68 : if( gpu ) {
392 : dline += " USEGPU";
393 : }
394 : } else {
395 : dline += " COMPONENTS";
396 : }
397 : // Parse origin
398 : std::string ostr;
399 147 : parse("ORIGIN",ostr);
400 : std::string num;
401 71 : if( ostr.length()>0 ) {
402 : // Parse atoms
403 24 : std::vector<std::string> afstr;
404 24 : MultiColvarShortcuts::parseAtomList("ATOMS",afstr,this);
405 17124 : for(unsigned i=0; i<afstr.size(); ++i) {
406 17100 : Tools::convert( i+1, num );
407 68400 : dline += " ATOMS" + num + "=" + ostr + "," + afstr[i];
408 : }
409 24 : } else {
410 47 : std::vector<std::string> grp;
411 47 : MultiColvarShortcuts::parseAtomList("GROUP",grp,this);
412 47 : std::vector<std::string> grpa;
413 94 : MultiColvarShortcuts::parseAtomList("GROUPA",grpa,this);
414 47 : if( grp.size()>0 ) {
415 3 : if( grpa.size()>0 ) {
416 0 : error("should not be using GROUPA in tandem with GROUP");
417 : }
418 : unsigned n=0;
419 203 : for(unsigned i=1; i<grp.size(); ++i) {
420 10103 : for(unsigned j=0; j<i; ++j) {
421 9903 : Tools::convert( n+1, num );
422 : n++;
423 39612 : dline += " ATOMS" + num + "=" + grp[i] + "," + grp[j];
424 : }
425 : }
426 44 : } else if( grpa.size()>0 ) {
427 1 : std::vector<std::string> grpb;
428 2 : MultiColvarShortcuts::parseAtomList("GROUPB",grpb,this);
429 1 : if( grpb.size()==0 ) {
430 0 : error("found GROUPA but no corresponding GROUPB");
431 : }
432 4 : for(unsigned i=0; i<grpa.size(); ++i) {
433 294 : for(unsigned j=0; j<grpb.size(); ++j) {
434 291 : Tools::convert( i*grpb.size() + j + 1, num );
435 1164 : dline += " ATOMS" + num + "=" + grpa[i] + "," + grpb[j];
436 : }
437 : }
438 1 : } else {
439 43 : std::string grpstr = getShortcutLabel() + "_grp: GROUP ATOMS=";
440 43 : for(unsigned i=1;; ++i) {
441 : std::string atstring;
442 1922 : parseNumbered("ATOMS",i,atstring);
443 961 : if( atstring.length()==0 ) {
444 : break;
445 : }
446 : std::string locstr;
447 1836 : parseNumbered("LOCATION",i,locstr);
448 918 : if( locstr.length()==0 ) {
449 301 : Tools::convert( i, num );
450 602 : readInputLine( getShortcutLabel() + "_vatom" + num + ": CENTER ATOMS=" + atstring );
451 301 : if( i==1 ) {
452 72 : grpstr += getShortcutLabel() + "_vatom" + num;
453 : } else {
454 795 : grpstr += "," + getShortcutLabel() + "_vatom" + num;
455 : }
456 : } else {
457 617 : if( i==1 ) {
458 : grpstr += locstr;
459 : } else {
460 610 : grpstr += "," + locstr;
461 : }
462 : }
463 918 : Tools::convert( i, num );
464 1836 : dline += " ATOMS" + num + "=" + atstring;
465 918 : }
466 43 : readInputLine( grpstr );
467 : }
468 47 : }
469 71 : readInputLine( dline );
470 : // Add shortcuts to label
471 71 : if( getName()=="DISTANCES" ) {
472 131 : MultiColvarShortcuts::expandFunctions( getShortcutLabel(), getShortcutLabel(), "", this );
473 3 : } else if( getName()=="XDISTANCES" ) {
474 2 : MultiColvarShortcuts::expandFunctions( getShortcutLabel(), getShortcutLabel() + ".x", "", this );
475 2 : } else if( getName()=="YDISTANCES" ) {
476 2 : MultiColvarShortcuts::expandFunctions( getShortcutLabel(), getShortcutLabel() + ".y", "", this );
477 1 : } else if( getName()=="ZDISTANCES" ) {
478 7 : MultiColvarShortcuts::expandFunctions( getShortcutLabel(), getShortcutLabel() + ".z", "", this );
479 : }
480 1604 : }
481 :
482 : }
483 : }
|
