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 "Between.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 BETWEEN
35 : /*
36 : Use a switching function to determine how many of the input variables are within a certain range.
37 :
38 : This action takes one argument, $s$ and evaluates the following function:
39 :
40 : $$
41 : w(s) = \int_a^b K\left( \frac{x - s}{\sigma} \right) \textrm{d}x
42 : $$
43 :
44 : In this equation $K$ is a symmetric function that must integrate to one and that is often
45 : called a [kernel function](https://en.wikipedia.org/wiki/Kernel_(statistics)) and $\sigma$ is a smearing parameter.
46 : The above function can be used to evaluate whether $s$ is between $a$ and $b$. The advantage of using the function above is that
47 : the resulting quantity has continuous derivatives. It can thus be used when calculating a collective variable upon which a simulation
48 : bias will be applied.
49 :
50 : The following example, shows how we can apply the function above on the instantaneous value of the distance between atom 1 and 2.
51 : The BETWEEN action here is used to determine whether the input distance is between 0.1 and 0.2 nm.
52 :
53 : ```plumed
54 : d: DISTANCE ATOMS=1,2
55 : b: BETWEEN ARG=d SWITCH={GAUSSIAN LOWER=0.1 UPPER=0.2 SMEAR=0.5}
56 : ```
57 :
58 : ## The Kernel function
59 :
60 : The $\sigma$ values in the expressions above is calculated from the parameters $a$, $b$ and $s$ that are provided using the
61 : `LOWER`, `UPPER` and `SMEAR` parameters respectively using the following function:
62 :
63 : $$
64 : \sigma = s(b-a)
65 : $$
66 :
67 : Also note that the Kernel function $K$ that is used in the first example is a Gaussian. The actual integral that is evaluated is thus:
68 :
69 : $$
70 : w(s) = \frac{1}{\sqrt{2\pi}\sigma} \int_a^b \exp\left( -\frac{(x-s)^2}{2\sigma^2}\right) \textrm{d}x
71 : $$
72 :
73 : The Gaussian kernel in this expression can be replaced by a triangular Kernel by changing the input to:
74 :
75 : ```plumed
76 : d: DISTANCE ATOMS=1,2
77 : b: BETWEEN ARG=d SWITCH={TRIANGULAR LOWER=0.1 UPPER=0.2 SMEAR=0.5}
78 : ```
79 :
80 : With this input the integral is then evaluated using:
81 :
82 : $$
83 : w(s) = \frac{1}{2\sigma} \int_a^b 1 - H\left( \frac{x-s}{\sigma} \right) \textrm{d}x
84 : $$
85 :
86 : where:
87 :
88 : $$
89 : H(x) = \begin{cases}
90 : x & \textrm{if} \quad x<1 \\
91 : 0 & \textrm{otherwise}
92 : \end{cases}
93 : $$
94 :
95 : ## Non rank zero arguments
96 :
97 : Instead of passing a single scalar in the input to the `BETWEEN` action you can pass a single vector as shown here:
98 :
99 : ```plumed
100 : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6 ATOMS4=7,8
101 : b: BETWEEN ARG=d SWITCH={GAUSSIAN LOWER=0.1 UPPER=0.2 SMEAR=0.5}
102 : ```
103 :
104 : The input to the `BETWEEN` action here is a vector with four elements. The output from the action `b` is similarly
105 : a vector with four elements. In calculating the elements of this vector PLUMED applies the function described in the previous
106 : 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 between
107 : 0.1 and 0.2 nm, the second element tells you if the distance between atoms 3 and 4 is between 0.1 and 0.2 nm and so on.
108 :
109 : You can use the commands in the above example input to evaluate the number of distances that are within the range of interest as follows:
110 :
111 : ```plumed
112 : d: DISTANCE ATOMS1=1,2 ATOMS2=3,4 ATOMS3=5,6 ATOMS4=7,8
113 : b: BETWEEN ARG=d SWITCH={GAUSSIAN LOWER=0.1 UPPER=0.2 SMEAR=0.5}
114 : s: SUM ARG=b PERIODIC=NO
115 : PRINT ARG=s FILE=colvar
116 : ```
117 :
118 : The final scalar that is output here is evaluated using the following summation:
119 :
120 : $$
121 : s = \sum_i \int_a^b \left( \frac{x - d_i}{\sigma} \right) \textrm{d}x
122 : $$
123 :
124 : 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
125 : between 0.1 and 0.2.
126 :
127 : Notice that you can do something similar with a matrix as input as shown below:
128 :
129 : ```plumed
130 : d: DISTANCE_MATRIX GROUPA=1-10 GROUPB=11-20
131 : b: BETWEEN ARG=d SWITCH={GAUSSIAN LOWER=0.1 UPPER=0.2 SMEAR=0.5}
132 : s: SUM ARG=b PERIODIC=NO
133 : PRINT ARG=s FILE=colvar
134 : ```
135 :
136 : 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
137 : tells you how many of these distances are between 0.1 and 0.2 nm.
138 :
139 : */
140 : //+ENDPLUMEDOC
141 :
142 : //+PLUMEDOC FUNCTION BETWEEN_VECTOR
143 : /*
144 : Use a switching function to determine how many of the input components are within a certain range
145 :
146 : \par Examples
147 :
148 : */
149 : //+ENDPLUMEDOC
150 :
151 : //+PLUMEDOC COLVAR BETWEEN_MATRIX
152 : /*
153 : Transform all the elements of a matrix using a switching function that is oen when the input value is within a particular range
154 :
155 : \par Examples
156 :
157 : */
158 : //+ENDPLUMEDOC
159 :
160 : typedef FunctionShortcut<Between> BetweenShortcut;
161 : PLUMED_REGISTER_ACTION(BetweenShortcut,"BETWEEN")
162 : typedef FunctionOfScalar<Between> ScalarBetween;
163 : PLUMED_REGISTER_ACTION(ScalarBetween,"BETWEEN_SCALAR")
164 : typedef FunctionOfVector<Between> VectorBetween;
165 : PLUMED_REGISTER_ACTION(VectorBetween,"BETWEEN_VECTOR")
166 : typedef FunctionOfMatrix<Between> MatrixBetween;
167 : PLUMED_REGISTER_ACTION(MatrixBetween,"BETWEEN_MATRIX")
168 :
169 215 : void Between::registerKeywords(Keywords& keys) {
170 215 : keys.add("compulsory","LOWER","the lower boundary for this particular bin");
171 215 : keys.add("compulsory","UPPER","the upper boundary for this particular bin");
172 215 : keys.add("compulsory","SMEAR","0.5","the ammount to smear the Gaussian for each value in the distribution");
173 215 : keys.add("optional","SWITCH","This keyword is used if you want to employ an alternative to the continuous function defined above. "
174 : "The following provides information on the \\ref histogrambead that are available. "
175 : "When this keyword is present you no longer need the LOWER, UPPER, SMEAR and KERNEL keywords.");
176 430 : keys.setValueDescription("scalar/vector/matrix","a function that is one if the input falls within a particular range and zero otherwise");
177 215 : }
178 :
179 53 : void Between::read( ActionWithArguments* action ) {
180 53 : if( action->getNumberOfArguments()!=1 ) {
181 0 : action->error("should only be one argument to between actions");
182 : }
183 :
184 : std::string str_min, str_max, tstr_min, tstr_max;
185 53 : bool isPeriodic = action->getPntrToArgument(0)->isPeriodic();
186 53 : if( isPeriodic ) {
187 1 : action->getPntrToArgument(0)->getDomain( str_min, str_max );
188 : }
189 :
190 : std::string hinput;
191 106 : action->parse("SWITCH",hinput);
192 53 : if(hinput.length()==0) {
193 : std::string low, up, sme;
194 5 : action->parse("LOWER",low);
195 5 : action->parse("UPPER",up);
196 5 : action->parse("SMEAR",sme);
197 10 : hinput = "GAUSSIAN LOWER=" + low + " UPPER=" + up + " SMEAR=" + sme;
198 : }
199 : std::string errors;
200 53 : hist.set( hinput, errors );
201 53 : if( errors.size()!=0 ) {
202 0 : action->error( errors );
203 : }
204 53 : action->log.printf(" %s \n", hist.description().c_str() );
205 :
206 53 : if( !isPeriodic ) {
207 : hist.isNotPeriodic();
208 : } else {
209 : double min;
210 1 : Tools::convert( str_min, min );
211 : double max;
212 1 : Tools::convert( str_max, max );
213 1 : hist.isPeriodic( min, max );
214 : }
215 53 : }
216 :
217 49226 : void Between::calc( const ActionWithArguments* action, const std::vector<double>& args, std::vector<double>& vals, Matrix<double>& derivatives ) const {
218 : plumed_dbg_assert( args.size()==1 );
219 49226 : vals[0] = hist.calculate( args[0], derivatives(0,0) );
220 49226 : }
221 :
222 : }
223 : }
224 :
225 :
|
