Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2014-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/ActionRegister.h"
23 : #include "tools/HistogramBead.h"
24 : #include "MultiColvarFilter.h"
25 :
26 : //+PLUMEDOC MTRANSFORMS MTRANSFORM_BETWEEN
27 : /*
28 : This action can be used to transform the colvar values calculated by a MultiColvar using a histogram bead
29 :
30 : In this action each colvar, \f$s_i\f$, calculated by MultiColvar is transformed by a \ref histogrambead function that
31 : is equal to one if the colvar is within a certain range and which is equal to zero otherwise. In other words, we
32 : compute:
33 :
34 : \f[
35 : f_i = \int_a^b K\left( \frac{s-s_i}{w} \right)
36 : \f]
37 :
38 : where \f$a, b\f$ and \f$w\f$ are parameters.
39 :
40 : It is important to understand the distinction between what is done here and what is done by \ref MFILTER_BETWEEN.
41 : In \ref MFILTER_BETWEEN a weight, \f$w_i\f$ for the colvar is calculated using the \ref histogrambead. If one calculates the
42 : MEAN for \ref MFILTER_BETWEEN one is thus calculating:
43 :
44 : \f[
45 : \mu = \frac{ \sum_i f_i s_i }{ \sum_i f_i}
46 : \f]
47 :
48 : In this action by contrast the colvar is being transformed by the \ref histogrambead. If one thus calculates a MEAN for
49 : this action one computes:
50 :
51 : \f[
52 : \mu = \frac{ \sum_{i=1}^N f_i }{ N }
53 : \f]
54 :
55 : In other words, you are calculating the mean for the transformed colvar.
56 :
57 : \par Examples
58 :
59 : The following input gives an example of how a \ref MTRANSFORM_BETWEEN action can be used to duplicate
60 : functionality that is elsewhere in PLUMED.
61 :
62 : \plumedfile
63 : DISTANCES ...
64 : GROUPA=1-10 GROUPB=11-20
65 : LABEL=d1
66 : ... DISTANCES
67 : MTRANSFORM_BETWEEN DATA=d1 LOWER=1.0 UPPER=2.0 SMEAR=0.5
68 : \endplumedfile
69 :
70 : In this case you can achieve the same result by using:
71 :
72 : \plumedfile
73 : DISTANCES ...
74 : GROUPA=1-10 GROUPB=11-20
75 : BETWEEN={GAUSSIAN LOWER=1.0 UPPER=2.0}
76 : ... DISTANCES
77 : \endplumedfile
78 : (see \ref DISTANCES)
79 :
80 : The advantage of \ref MTRANSFORM_BETWEEN comes, however, if you want to use transformed colvars as input
81 : for \ref MULTICOLVARDENS
82 :
83 : */
84 : //+ENDPLUMEDOC
85 :
86 : //+PLUMEDOC MFILTERS MFILTER_BETWEEN
87 : /*
88 : This action can be used to filter the colvar values calculated by a multicolvar so that one can compute the mean and so on for only those multicolvars within a certain range.
89 :
90 : This action can be used to create a dynamic group of atom based on the value of a multicolvar.
91 : In this action a multicolvar is within the dynamic group if its value lies in a particular range.
92 : In actuality a weight, \f$w_i\f$ is ascribed to each colvar, \f$s_i\f$ calculated by a multicolvar
93 : and this weight measures the degree to which a colvar is a member of the group. This weight is
94 : calculated using a \ref histogrambead so it is given by:
95 :
96 : \f[
97 : w_i = \int_a^b K\left( \frac{s - s_i}{w} \right)
98 : \f]
99 :
100 : where \f$a, b\f$ and \f$w\f$ are parameters. If one calculates a function of the set of multicolvars
101 : these weights are included in the calculation. As such if one calculates the MEAN, \f$\mu\f$ of a filtered
102 : multicolvar what is computed is the following:
103 :
104 : \f[
105 : \mu = \frac{ \sum_i w_i s_i }{ \sum_i w_i}
106 : \f]
107 :
108 : One is thus calculating the mean for those colvars that are within the range of interest.
109 :
110 : \par Examples
111 :
112 : The example shown below calculates the mean for those distances that are between 0 and 3 nm in length
113 :
114 : \plumedfile
115 : DISTANCES GROUPA=1 GROUPB=2-50 MEAN LABEL=d1
116 : MFILTER_BETWEEN DATA=d1 LOWER=0 UPPER=3.0 SMEAR=0.0001 MEAN LABEL=d4
117 : \endplumedfile
118 :
119 : More complicated things can be done by using the label of a filter as input to a new multicolvar as shown
120 : in the example below. Here the coordination numbers of all atoms are computed. The atoms with a coordination
121 : number between 4 and 6 are then identified using the filter. This reduced list of atoms is then used as input
122 : to a second coordination number calculation. This second coordination number thus measures the number of atoms
123 : 4-6 coordinated atoms each of the 4-6 coordination atoms is bound to.
124 :
125 : \plumedfile
126 : c1: COORDINATIONNUMBER SPECIES=1-150 SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0}
127 : cf: MFILTER_BETWEEN DATA=c1 LOWER=4 UPPER=6 SMEAR=0.5 LOWMEM
128 : c2: COORDINATIONNUMBER SPECIES=cf SWITCH={EXP D_0=4.0 R_0=0.5 D_MAX=6.0} MORE_THAN={RATIONAL D_0=2.0 R_0=0.1}
129 : \endplumedfile
130 :
131 : */
132 : //+ENDPLUMEDOC
133 :
134 : namespace PLMD {
135 : namespace multicolvar {
136 :
137 : class FilterBetween : public MultiColvarFilter {
138 : private:
139 : HistogramBead hb;
140 : public:
141 : static void registerKeywords( Keywords& keys );
142 : explicit FilterBetween(const ActionOptions& ao);
143 : double applyFilter( const double& val, double& df ) const override;
144 : };
145 :
146 13791 : PLUMED_REGISTER_ACTION(FilterBetween,"MFILTER_BETWEEN")
147 13785 : PLUMED_REGISTER_ACTION(FilterBetween,"MTRANSFORM_BETWEEN")
148 :
149 11 : void FilterBetween::registerKeywords( Keywords& keys ) {
150 11 : MultiColvarFilter::registerKeywords( keys );
151 22 : keys.add("compulsory","LOWER","the lower boundary for the range of interest");
152 22 : keys.add("compulsory","UPPER","the upper boundary for the range of interest");
153 22 : keys.add("compulsory","SMEAR","0.5","the amount by which to smear the value for kernel density estimation");
154 22 : keys.add("optional","BEAD","This keywords is used if you want to employ an alternative to the function defined above. "
155 : "The following provides information on the \\ref histogrambead that are available. "
156 : "When this keyword is present you no longer need the LOWER, UPPER and SMEAR keywords.");
157 11 : }
158 :
159 3 : FilterBetween::FilterBetween(const ActionOptions& ao):
160 : Action(ao),
161 3 : MultiColvarFilter(ao) {
162 : // Read in the switching function
163 : std::string sw, errors;
164 6 : parse("BEAD",sw);
165 3 : if( getPntrToMultiColvar()->isPeriodic() ) {
166 : std::string min, max;
167 0 : getPntrToMultiColvar()->retrieveDomain( min, max );
168 : double mlow, mhigh;
169 0 : Tools::convert( min,mlow );
170 0 : Tools::convert( max, mhigh);
171 0 : hb.isPeriodic( mlow, mhigh );
172 : } else {
173 : hb.isNotPeriodic();
174 : }
175 :
176 3 : if(sw.length()>0) {
177 0 : hb.set(sw,errors);
178 0 : if( errors.length()!=0 ) {
179 0 : error("problem reading BEAD keyword : " + errors );
180 : }
181 : } else {
182 : double l, u, s;
183 : std::string ll, uu, ss;
184 3 : parse("LOWER",l);
185 3 : parse("UPPER",u);
186 3 : parse("SMEAR",s);
187 3 : Tools::convert(l,ll);
188 3 : Tools::convert(u,uu);
189 3 : Tools::convert(s,ss);
190 6 : sw="GAUSSIAN LOWER=" + ll + " UPPER=" + uu + " SMEAR=" + ss;
191 3 : hb.set(sw,errors);
192 3 : plumed_massert(errors.length()==0,"problems with bead" + errors);
193 : }
194 3 : log.printf(" filtering colvar values and focussing only on those values in range %s\n",( hb.description() ).c_str() );
195 :
196 3 : checkRead();
197 3 : }
198 :
199 735 : double FilterBetween::applyFilter( const double& val, double& df ) const {
200 735 : double f = hb.calculate( val, df );
201 735 : return f;
202 : }
203 :
204 : }
205 : }
|
