Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-2017 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 "LessThan.h"
23 : #include "FunctionShortcut.h"
24 : #include "FunctionOfScalar.h"
25 : #include "FunctionOfVector.h"
26 : #include "FunctionOfMatrix.h"
27 : #include "core/ActionRegister.h"
28 :
29 : #include <cmath>
30 :
31 : namespace PLMD {
32 : namespace function {
33 :
34 : //+PLUMEDOC FUNCTION LESS_THAN
35 : /*
36 : Use a switching function to determine how many of the input variables are less than a certain cutoff.
37 :
38 : This action takes one argument, $r$ and evaluates the following function:
39 :
40 : $$
41 : w(r) = s(r)
42 : $$
43 :
44 : The function $s(r)$ here is a switching function so the output value $w$ is a number between 0 and 1.
45 : Switching functions are typically used to determine if an input value is less than some cutoff. The value
46 : of $w$ the switches smoothly from one to zero as the input value $r$ crosses the threshold of interest.
47 :
48 : The following example, shows how we can apply a switching function on the instantaneous value of the
49 : distance between atom 1 and atom 2.
50 :
51 : ```plumed
52 : d: DISTANCE ATOMS=1,2
53 : s: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.2}
54 : ```
55 :
56 : The output here is close to one when the distance between atoms 1 and 2 is less than 0.2 nm close to zero
57 : for values that are greater than 0.2.
58 :
59 : ## Non rank zero arguments
60 :
61 : Instead of passing a single scalar in the input to the `LESS_THAN` action you can pass a single vector as shown here:
62 :
63 : ```plumed
64 : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6 ATOMS4=7,8
65 : b: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.2}
66 : ```
67 :
68 : The input to the `LESS_THAN` action here is a vector with four elements. The output from the action `b` is similarly
69 : a vector with four elements. In calculating the elements of this vector PLUMED applies the function described in the previous
70 : section on each of the distances in turn. The first element of `b` thus tells you if the distance between atoms 1 and 2 is less than
71 : 0.2 nm, the second element tells you if the distance between atoms 3 and 4 is less than 0.2 nm and so on.
72 :
73 : You can use the commands in the above example input to evaluate the number of distances that are less than a cutoff as follows:
74 :
75 : ```plumed
76 : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6 ATOMS4=7,8
77 : b: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.2}
78 : s: SUM ARG=b PERIODIC=NO
79 : PRINT ARG=s FILE=colvar
80 : ```
81 :
82 : The final scalar that is output here is evaluated using the following summation:
83 :
84 : $$
85 : s = \sum_i s(d_i)
86 : $$
87 :
88 : where the sum over $i$ here runs over the four distances in the above expression. This scalar tells you the number of distances that are
89 : less than 0.2 nm.
90 :
91 : Notice that you can do something similar with a matrix as input as shown below:
92 :
93 : ```plumed
94 : d: DISTANCE_MATRIX GROUPA=1-10 GROUPB=11-20
95 : b: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.2}
96 : s: SUM ARG=b PERIODIC=NO
97 : PRINT ARG=s FILE=colvar
98 : ```
99 :
100 : This input tells PLUMED to calculate the 100 distances between the atoms in the two input groups. The final value that is printed to the colvar file then
101 : tells you how many of these distances are less than 0.2 nm.
102 :
103 : ## Switching functions types
104 :
105 : PLUMED allows you to use a range of different switching function types. Most of these
106 : require you to provide at least one input parameter $r_0$. The switching function is then defined so that for
107 : $r \le d_0 \quad s(r)=1.0$ while for $r > d_0$ the function decays smoothly to 0. By changing the switching
108 : function you are using you change the decay that occurs in the switching function for $r$ values that are greater
109 : that $d_0$.
110 :
111 : You can specify the value at which the the switching function goes to zero by using the D_MAX keyword.
112 : PLUMED will ensure that the swtiching function is definitely zero for for input values that are greater than
113 : or equal to D_MAX by by stretching and shifting the function using the following transformations.
114 :
115 : $
116 : s(r) = \frac{s'(r)-s'(d_{max})}{s'(0)-s'(d_{max})}
117 : $
118 :
119 : The various expressions that are used for $s'(r)$ in the above expression are described in the following sections.
120 : Scaling the switching function so that it is exactly zero at $d_{max}$ using the transformation described above
121 : has been the default behaviour within PLUMED since version 2.2.
122 : If you wish to use a D_MAX parameter and the $s'(r)$ functions described in the following sections directly you can use the
123 : keyword `NOSTRETCH`. However, the advangage of performing the stretching and scaling in this way is that the function has
124 : no discontinuities. If you use the `NOSTRETCH` option the switching function may have a discontinuity at $d_{max}$.
125 :
126 : ### Rational switching function
127 :
128 : The rational switching function is the most commonly used of the switching functions. The implementation of this function within PLUMED has been carefully
129 : optimised so calculations using this switching function should be fast.
130 :
131 : For this function $s'(r)$ has the following functional form:
132 :
133 : $$
134 : s'(r) = \begin{cases}
135 : 1 & \textrm{if} \quad r<d_0 \\
136 : \frac{ 1 - \left(\frac{ r - d_0 }{ r_0 }\right)^{n} }{ 1 - \left(\frac{ r - d_0 }{ r_0 }\right)^{m} } & \textrm{if} \quad d_0 < r < d_{max} \\
137 : 0 & \textrm{otherwise}
138 : $$
139 :
140 : The following example illustrates how you can use the rational switching function is used in practice:
141 :
142 : ```plumed
143 : d: DISTANCE ATOMS=1,2
144 : s: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.5 D_0=0.1 NN=8 MM=16 D_MAX=1.0}
145 : ```
146 :
147 : As was discussed in earlier parts of this page, you can also specify a rational switching function using the following input:
148 :
149 : ```plumed
150 : d: DISTANCE ATOMS=1,2
151 : s: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.2}
152 : ```
153 :
154 : In this case, the paramers $d_0$, $n$ and $m$ are set to their default values of 0, 6 and 12 respectively. Meanwhile, if D_MAX is unset the stretching and scaling that
155 : was described in the previous section is not performed.
156 :
157 : Notice that you can choose to set a subset of the switching function parameters and to use the default values for the unset parameters by using an input like this one:
158 :
159 : ```plumed
160 : d: DISTANCE ATOMS=1,2
161 : s: LESS_THAN ARG=d SWITCH={RATIONAL R_0=0.2 NN=8}
162 : ```
163 :
164 : The input above sets $n=8$ and $m=2n=16$. Meanwhile, $d_0$ is set to its default value of 0.
165 :
166 : Lastly, note that you can also specify the parameters for the rational switching function without using the SWITCH keyword as indicated below:
167 :
168 : ```plumed
169 : d: DISTANCE ATOMS=1,2
170 : s: LESS_THAN ARG=d D_0=0.1 R_0=0.2 NN=6 MM=12
171 : ```
172 :
173 : The newer syntax using the SWITCH keyword is better than this option as there is no way to set the D_MAX parameter this way and because if you use this syntax you are forced to
174 : use the RATIONAL switching function type.
175 :
176 : ### Expoential switching function
177 :
178 : The way that the exponential switching function can be used is illustrated in the following example input:
179 :
180 : ```plumed
181 : d: DISTANCE ATOMS=1,2
182 : s: LESS_THAN ARG=d SWITCH={EXP D_0=0.1 R_0=0.2 D_MAX=1.0}
183 : ```
184 :
185 : The $s'(r)$ that is used in this input is given by the following expression:
186 :
187 : $$
188 : s'(r) = \begin{cases}
189 : 1 & \textrm{if} \quad r<d_0 \\
190 : \exp\left(-\frac{ r - d_0 }{ r_0 }\right) & \textrm{if} \quad d_0 < r < d_{max} \\
191 : 0 & \textrm{otherwise}
192 : $$
193 :
194 : You can also specify that an exponential switching function is to be used by using the following input:
195 :
196 : ```plumed
197 : d: DISTANCE ATOMS=1,2
198 : s: LESS_THAN ARG=d SWITCH={EXP R_0=0.2}
199 : ```
200 :
201 : For this input $d_0$ is set to its default value of 0. Furthermore, as D_MAX is unset the stretching and scaling that
202 : was described in the previous section is not performed.
203 :
204 : ### Gaussian switching function
205 :
206 : The way that the gaussian switching function can be used is illustrated in the following example input:
207 :
208 : ```plumed
209 : d: DISTANCE ATOMS=1,2
210 : s: LESS_THAN ARG=d SWITCH={GAUSSIAN D_0=0.1 R_0=0.2 D_MAX=1.0}
211 : ```
212 :
213 : The $s'(r)$ that is used in this input is given by the following expression:
214 :
215 : $$
216 : s'(r) = \begin{cases}
217 : 1 & \textrm{if} \quad r<d_0 \\
218 : \exp\left(-\frac{ (r - d_0)^2 }{ 2r_0^2 }\right) & \textrm{if} \quad d_0 < r < d_{max} \\
219 : 0 & \textrm{otherwise}
220 : $$
221 :
222 : You can also specify that an gaussian switching function is to be used by using the following input:
223 :
224 : ```plumed
225 : d: DISTANCE ATOMS=1,2
226 : s: LESS_THAN ARG=d SWITCH={GAUSSIAN R_0=0.2}
227 : ```
228 :
229 : For this input $d_0$ is set to its default value of 0. Furthermore, as D_MAX is unset the stretching and scaling that
230 : was described in the previous section is not performed.
231 :
232 : ### Sketch-map switching function
233 :
234 : You can specify that the sketch-map switching function is to be used by using the following input:
235 :
236 : ```plumed
237 : d: DISTANCE ATOMS=1,2
238 : s: LESS_THAN ARG=d SWITCH={SMAP D_0=1 R_0=4 A=3 B=2 D_MAX=10}
239 : ```
240 :
241 : The $s'(r)$ that is used in this input is given by the following expression:
242 :
243 : $$
244 : s'(r) = \begin{cases}
245 : 1 & \textrm{if} \quad r<d_0 \\
246 : \left[ 1 + ( 2^{a/b} -1 )\left( \frac{r-d_0}{r_0} \right)^a \right]^{-b/a} & \textrm{if} \quad d_0 < r < d_{max} \\
247 : 0 & \textrm{otherwise}
248 : $$
249 :
250 : This type of swtiching function would be used with PLUMED's implementation of [SKETCHMAP](SKETCHMAP.md). If you are
251 : performing a sketch-map projection an input like the one below would be more common:
252 :
253 : ```plumed
254 : d: DISTANCE ATOMS=1,2
255 : s: LESS_THAN ARG=d SWITCH={SMAP R_0=4 A=3 B=2}
256 : ```
257 :
258 : With this input $d_0$ is set to its default value of 0. Furthermore, as D_MAX is unset the stretching and scaling that
259 : was described in the previous section is not performed.
260 :
261 : ### Q switching function
262 :
263 : You can specify that the Q switching function is to be used by using the following input:
264 :
265 : ```plumed
266 : d: DISTANCE ATOMS=1,2
267 : s: LESS_THAN ARG=d SWITCH={Q D_0=0 REF=0.498366 BETA=50.0 LAMBDA=1.8 R_0=0.01 D_MAX=20}
268 : ```
269 :
270 : The $s'(r)$ that is used in this input is given by the following expression:
271 :
272 : $$
273 : s'(r) = \begin{cases}
274 : 1 & \textrm{if} \quad r<d_0 \\
275 : s(r) = \frac{1}{1 + \exp(\beta(r_{ij} - \lambda r_{ij}^0))} & \textrm{if} \quad d_0 < r < d_{max} \\
276 : 0 & \textrm{otherwise}
277 : $$
278 :
279 : The value of $r_{ij}^0$ is specified in the above input by the `REF` keyword. Notice that you can also specify this type of switching function using
280 : the following input:
281 :
282 : ```plumed
283 : d: DISTANCE ATOMS=1,2
284 : s: LESS_THAN ARG=d SWITCH={Q REF=0.498366 R_0=0.01 }
285 : ```
286 :
287 : With this input $_0$, $\lambda$ and $\beta$ are set equal to their default values of 0, 1.8 and 50 nm$^{-1}$ respectively. Furthermore, as D_MAX is unset the stretching and scaling that
288 : was described in the previous section is not performed. These default values for $\lambda$ and $\beta$ are suitable if you are simulating all atom models. However, if you are using
289 : a coarse grainded models values of $\lambda=1.5$ and $\beta=50 nm^-1$ are more appropriate.
290 :
291 : ### Cubic switching function
292 :
293 : You can specify that the cubic swtiching function is to be used by using the following input:
294 :
295 : ```plumed
296 : d: DISTANCE ATOMS=1,2
297 : s: LESS_THAN ARG=d SWITCH={CUBIC D_0=0.1 D_MAX=0.5}
298 : ```
299 :
300 : If this type of expression is used then the $s(r)$ is calculated as:
301 :
302 : $$
303 : s(r) = (y-1)^2(1+2y) \qquad \textrm{where} \quad y = \frac{d_0 - d_{max}}{d_0-d_{max}}$
304 : $$
305 :
306 : No stretching is required for this type of switching function as its functional form ensures that it is zero at $d_{max}$.
307 :
308 : ### Tanh switching function
309 :
310 : You can specify that the cubic swtiching function is to be used by using the following input:
311 :
312 : ```plumed
313 : d: DISTANCE ATOMS=1,2
314 : s: LESS_THAN ARG=d SWITCH={TANH D_0=0.1 R_0=0.1 D_MAX=1.0}
315 : ```
316 :
317 : The $s'(r)$ that is used in this input is given by the following expression:
318 :
319 : $$
320 : s'(r) = \begin{cases}
321 : 1 & \textrm{if} \quad r<d_0 \\
322 : 1 - \tanh\left( \frac{ r - d_0 }{ r_0 } \right) & \textrm{if} \quad d_0 < r < d_{max} \\
323 : 0 & \textrm{otherwise}
324 : $$
325 :
326 : You can also specify that a tanh switching function is to be used by using the following input:
327 :
328 : ```plumed
329 : d: DISTANCE ATOMS=1,2
330 : s: LESS_THAN ARG=d SWITCH={TANH R_0=0.2}
331 : ```
332 :
333 : For this input $d_0$ is set to its default value of 0. Furthermore, as D_MAX is unset the stretching and scaling that
334 : was described in the previous section is not performed.
335 :
336 : ### Cosinus switching function
337 :
338 : You can specify that the cosinum function is to be used by using the following input:
339 :
340 : ```plumed
341 : d: DISTANCE ATOMS=1,2
342 : s: LESS_THAN ARG=d SWITCH={COSINUS D_0=0.1 R_0=0.1}
343 : ```
344 :
345 : The $s'(r)$ that is used in this input is given by the following expression:
346 :
347 : $$
348 : s(r) =\left\{\begin{array}{ll}
349 : 1 & \mathrm{if } r \leq d_0 \\
350 : 0.5 \left( \cos ( \frac{ r - d_0 }{ r_0 } \pi ) + 1 \right) & \mathrm{if } d_0 < r\leq d_0 + r_0 \\
351 : 0
352 : $$
353 :
354 : No stretching is required for this type of switching function as its functional form ensures that it is zero at $d_0+r_0$
355 :
356 : ### Custom switching function
357 :
358 : If you want to use a switching function that is different to all of the functions listed above you can use the cusom option that is illustrated below:
359 :
360 : ```plumed
361 : d: DISTANCE ATOMS=1,2
362 : s: LESS_THAN ARG=d SWITCH={CUSTOM FUNC=1/(1+x^6) D_0=0.1 R_0=0.1 D_MAX=10}
363 : ```
364 :
365 : This option allows you to use the lepton library that is used to implement the [CUSTOM](CUSTOM.md) action to evaluate your switching function. With this
366 : option you specify the functional form for $s'(r)$ by providing a function of `x` to the `FUNC` keyword. The `x` that enters the switching function
367 : definition that you provide is given by:
368 :
369 : $$
370 : x = \frac{r - d_0}{r_0}
371 : $$
372 :
373 : Notice, that you can also express your switching function in terms of $x^2$ as illustrated in the input below:
374 :
375 : ```plumed
376 : d: DISTANCE ATOMS=1,2
377 : s: LESS_THAN ARG=d SWITCH={CUSTOM FUNC=1/(1+x2^3) D_0=0.1 R_0=0.1 D_MAX=10}
378 : ```
379 :
380 : which is equivalent to the earlier input with the CUSTOM switching function. However, using `x2` in place of `x` is often more computationally efficient
381 : as you can avoid an expensive square root operation. Notice, lastly, that just as there is a [MATHEVAL](MATHEVAL.md) action that is equivalent to [CUSTOM](CUSTOM.md),
382 : there is also a matheval equivalent that you can use here:
383 :
384 : ```plumed
385 : d: DISTANCE ATOMS=1,2
386 : s: LESS_THAN ARG=d SWITCH={MATHEVAL FUNC=1/(1+x2^3) R_0=0.1}
387 : ```
388 :
389 : For this input $d_0$ is set to its default value of 0. Furthermore, as D_MAX is unset the stretching and scaling that
390 : was described in the previous section is not performed.
391 :
392 : > [!CAUTION]
393 : > With the default implementation CUSTOM is slower than other functions
394 : > (e.g., it is slower than an equivalent RATIONAL function by approximately a factor 2).
395 : > You can find information on how to improve its performance in the documenation for [CUSTOM](CUSTOM.md)
396 :
397 : */
398 : //+ENDPLUMEDOC
399 :
400 : //+PLUMEDOC FUNCTION LESS_THAN_VECTOR
401 : /*
402 : Use a switching function to determine how many components of the vector are less than a certain cutoff.
403 :
404 : \par Examples
405 :
406 : */
407 : //+ENDPLUMEDOC
408 :
409 : //+PLUMEDOC COLVAR LESS_THAN_MATRIX
410 : /*
411 : Transform all the elements of a matrix using a switching function that is one when the input value is smaller than a threshold
412 :
413 : \par Examples
414 :
415 : */
416 : //+ENDPLUMEDOC
417 :
418 : typedef FunctionShortcut<LessThan> LessThanShortcut;
419 : PLUMED_REGISTER_ACTION(LessThanShortcut,"LESS_THAN")
420 : typedef FunctionOfScalar<LessThan> ScalarLessThan;
421 : PLUMED_REGISTER_ACTION(ScalarLessThan,"LESS_THAN_SCALAR")
422 : typedef FunctionOfVector<LessThan> VectorLessThan;
423 : PLUMED_REGISTER_ACTION(VectorLessThan,"LESS_THAN_VECTOR")
424 : typedef FunctionOfMatrix<LessThan> MatrixLessThan;
425 : PLUMED_REGISTER_ACTION(MatrixLessThan,"LESS_THAN_MATRIX")
426 :
427 300 : void LessThan::registerKeywords(Keywords& keys) {
428 300 : keys.add("compulsory","NN","6","The n parameter of the switching function ");
429 300 : keys.add("compulsory","MM","0","The m parameter of the switching function; 0 implies 2*NN");
430 300 : keys.add("compulsory","D_0","0.0","The d_0 parameter of the switching function");
431 300 : keys.add("compulsory","R_0","The r_0 parameter of the switching function");
432 300 : keys.add("optional","SWITCH","This keyword is used if you want to employ an alternative to the continuous swiching function defined above. "
433 : "The following provides information on the \\ref switchingfunction that are available. "
434 : "When this keyword is present you no longer need the NN, MM, D_0 and R_0 keywords.");
435 300 : keys.addFlag("SQUARED",false,"is the input quantity the square of the value that you would like to apply the switching function to");
436 600 : keys.setValueDescription("scalar/vector/matrix","a function that is one if the input is less than a threshold");
437 300 : }
438 :
439 73 : void LessThan::read( ActionWithArguments* action ) {
440 73 : if( action->getNumberOfArguments()!=1 ) {
441 0 : action->error("should only be one argument to less_than actions");
442 : }
443 73 : if( action->getPntrToArgument(0)->isPeriodic() ) {
444 0 : action->error("cannot use this function on periodic functions");
445 : }
446 :
447 :
448 : std::string sw,errors;
449 146 : action->parse("SWITCH",sw);
450 73 : if(sw.length()>0) {
451 73 : switchingFunction.set(sw,errors);
452 73 : if( errors.length()!=0 ) {
453 0 : action->error("problem reading SWITCH keyword : " + errors );
454 : }
455 : } else {
456 0 : int nn=6;
457 0 : int mm=0;
458 0 : double d0=0.0;
459 0 : double r0=0.0;
460 0 : action->parse("R_0",r0);
461 0 : if(r0<=0.0) {
462 0 : action->error("R_0 should be explicitly specified and positive");
463 : }
464 0 : action->parse("D_0",d0);
465 0 : action->parse("NN",nn);
466 0 : action->parse("MM",mm);
467 0 : switchingFunction.set(nn,mm,r0,d0);
468 : }
469 73 : action->log<<" using switching function with cutoff "<<switchingFunction.description()<<"\n";
470 73 : action->parseFlag("SQUARED",squared);
471 73 : if( squared ) {
472 0 : action->log<<" input quantity is square of quantity that switching function acts upon\n";
473 : }
474 73 : }
475 :
476 203171 : void LessThan::calc( const ActionWithArguments* action, const std::vector<double>& args, std::vector<double>& vals, Matrix<double>& derivatives ) const {
477 : plumed_dbg_assert( args.size()==1 );
478 203171 : if( squared ) {
479 0 : vals[0] = switchingFunction.calculateSqr( args[0], derivatives(0,0) );
480 : } else {
481 203171 : vals[0] = switchingFunction.calculate( args[0], derivatives(0,0) );
482 : }
483 203171 : derivatives(0,0) = args[0]*derivatives(0,0);
484 203171 : }
485 :
486 : }
487 : }
488 :
489 :
|
