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 : #include "MultiColvarShortcuts.h"
23 : #include "core/PlumedMain.h"
24 : #include "core/ActionSet.h"
25 : #include "core/Group.h"
26 :
27 : namespace PLMD {
28 : namespace multicolvar {
29 :
30 720 : void MultiColvarShortcuts::shortcutKeywords( Keywords& keys ) {
31 720 : keys.add("numbered","LESS_THAN","calculate the number of variables that are less than a certain target value. "
32 : "This quantity is calculated using \\f$\\sum_i \\sigma(s_i)\\f$, where \\f$\\sigma(s)\\f$ "
33 : "is a \\ref switchingfunction.");
34 1440 : keys.linkActionInDocs("LESS_THAN","LESS_THAN");
35 1440 : keys.reset_style("LESS_THAN","deprecated");
36 1440 : keys.addOutputComponent("lessthan","LESS_THAN","scalar","the number of colvars that have a value less than a threshold");
37 720 : keys.add("numbered","MORE_THAN","calculate the number of variables that are more than a certain target value. "
38 : "This quantity is calculated using \\f$\\sum_i 1 - \\sigma(s_i)\\f$, where \\f$\\sigma(s)\\f$ "
39 : "is a \\ref switchingfunction.");
40 1440 : keys.linkActionInDocs("MORE_THAN","MORE_THAN");
41 1440 : keys.reset_style("MORE_THAN","deprecated");
42 1440 : keys.addOutputComponent("morethan","MORE_THAN","scalar","the number of colvars that have a value more than a threshold");
43 720 : keys.add("optional","ALT_MIN","calculate the minimum value. "
44 : "To make this quantity continuous the minimum is calculated using "
45 : "\\f$ \\textrm{min} = -\\frac{1}{\\beta} \\log \\sum_i \\exp\\left( -\\beta s_i \\right) \\f$ "
46 : "The value of \\f$\\beta\\f$ in this function is specified using (BETA=\\f$\\beta\\f$).");
47 1440 : keys.reset_style("ALT_MIN","deprecated");
48 1440 : keys.addOutputComponent("altmin","ALT_MIN","scalar","the minimum value of the cv");
49 720 : keys.add("optional","MIN","calculate the minimum value. "
50 : "To make this quantity continuous the minimum is calculated using "
51 : "\\f$ \\textrm{min} = \\frac{\\beta}{ \\log \\sum_i \\exp\\left( \\frac{\\beta}{s_i} \\right) } \\f$ "
52 : "The value of \\f$\\beta\\f$ in this function is specified using (BETA=\\f$\\beta\\f$)");
53 1440 : keys.reset_style("MIN","deprecated");
54 1440 : keys.addOutputComponent("min","MIN","scalar","the minimum colvar");
55 720 : keys.add("optional","MAX","calculate the maximum value. "
56 : "To make this quantity continuous the maximum is calculated using "
57 : "\\f$ \\textrm{max} = \\beta \\log \\sum_i \\exp\\left( \\frac{s_i}{\\beta}\\right) \\f$ "
58 : "The value of \\f$\\beta\\f$ in this function is specified using (BETA=\\f$\\beta\\f$)");
59 1440 : keys.reset_style("MAX","deprecated");
60 1440 : keys.addOutputComponent("max","MAX","scalar","the maximum colvar");
61 720 : keys.add("numbered","BETWEEN","calculate the number of values that are within a certain range. "
62 : "These quantities are calculated using kernel density estimation as described on "
63 : "\\ref histogrambead.");
64 1440 : keys.linkActionInDocs("BETWEEN","BETWEEN");
65 1440 : keys.reset_style("BETWEEN","deprecated");
66 1440 : keys.addOutputComponent("between","BETWEEN","scalar","the number of colvars that have a value that lies in a particular interval");
67 720 : keys.addFlag("HIGHEST",false,"this flag allows you to recover the highest of these variables.");
68 1440 : keys.reset_style("HIGHEST","deprecated");
69 1440 : keys.addOutputComponent("highest","HIGHEST","scalar","the largest of the colvars");
70 720 : keys.add("optional","HISTOGRAM","calculate a discretized histogram of the distribution of values. "
71 : "This shortcut allows you to calculates NBIN quantites like BETWEEN.");
72 1440 : keys.reset_style("HISTOGRAM","deprecated");
73 720 : keys.addFlag("LOWEST",false,"this flag allows you to recover the lowest of these variables.");
74 1440 : keys.reset_style("LOWEST","deprecated");
75 1440 : keys.addOutputComponent("lowest","LOWEST","scalar","the smallest of the colvars");
76 720 : keys.addFlag("SUM",false,"calculate the sum of all the quantities.");
77 1440 : keys.reset_style("SUM","deprecated");
78 1440 : keys.addOutputComponent("sum","SUM","scalar","the sum of the colvars");
79 720 : keys.addFlag("MEAN",false,"calculate the mean of all the quantities.");
80 1440 : keys.reset_style("MEAN","deprecated");
81 1440 : keys.addOutputComponent("mean","MEAN","scalar","the mean of the colvars");
82 720 : keys.needsAction("SUM");
83 720 : keys.needsAction("MEAN");
84 720 : keys.needsAction("CUSTOM");
85 720 : keys.needsAction("HIGHEST");
86 720 : keys.needsAction("LOWEST");
87 720 : keys.needsAction("LESS_THAN");
88 720 : keys.needsAction("MORE_THAN");
89 720 : keys.needsAction("BETWEEN");
90 720 : }
91 :
92 119 : void MultiColvarShortcuts::expandFunctions( const std::string& labout, const std::string& argin, const std::string& weights, ActionShortcut* action ) {
93 : std::map<std::string,std::string> keymap;
94 119 : readShortcutKeywords( keymap, action );
95 119 : expandFunctions( labout, argin, weights, keymap, action );
96 119 : }
97 :
98 234 : void MultiColvarShortcuts::readShortcutKeywords( std::map<std::string,std::string>& keymap, ActionShortcut* action ) {
99 234 : Keywords keys;
100 234 : shortcutKeywords( keys );
101 234 : action->readShortcutKeywords( keys, keymap );
102 234 : }
103 :
104 102 : void MultiColvarShortcuts::parseAtomList( const std::string& key, std::vector<std::string>& atoms, ActionShortcut* action ) {
105 : std::vector<std::string> astr;
106 102 : action->parseVector(key,astr);
107 102 : if( astr.size()==0 ) {
108 : return ;
109 : }
110 36 : Tools::interpretRanges( astr );
111 1281 : for(unsigned i=0; i<astr.size(); ++i) {
112 1245 : Group* mygr=action->plumed.getActionSet().selectWithLabel<Group*>(astr[i]);
113 1245 : if( mygr ) {
114 5 : std::vector<std::string> grstr( mygr->getGroupAtoms() );
115 121 : for(unsigned j=0; j<grstr.size(); ++j) {
116 116 : atoms.push_back(grstr[j]);
117 : }
118 5 : } else {
119 1240 : Group* mygr2=action->plumed.getActionSet().selectWithLabel<Group*>(astr[i] + "_grp");
120 1240 : if( mygr2 ) {
121 10 : std::vector<std::string> grstr( mygr2->getGroupAtoms() );
122 16082 : for(unsigned j=0; j<grstr.size(); ++j) {
123 16072 : atoms.push_back(grstr[j]);
124 : }
125 10 : } else {
126 1230 : atoms.push_back(astr[i]);
127 : }
128 : }
129 : }
130 102 : }
131 :
132 243 : void MultiColvarShortcuts::expandFunctions( const std::string& labout, const std::string& argin, const std::string& weights,
133 : const std::map<std::string,std::string>& keymap, ActionShortcut* action ) {
134 243 : if( keymap.empty() ) {
135 : return;
136 : }
137 : // Parse LESS_THAN
138 294 : if( keymap.count("LESS_THAN") ) {
139 16 : std::string sum_arg = labout + "_lt", lt_string = keymap.find("LESS_THAN")->second;
140 16 : action->readInputLine( labout + "_lt: LESS_THAN ARG=" + argin + " SWITCH={" + lt_string + "}");
141 8 : if( weights.length()>0 ) {
142 1 : sum_arg = labout + "_wlt";
143 2 : action->readInputLine( labout + "_wlt: CUSTOM ARG=" + weights + "," + labout + "_lt FUNC=x*y PERIODIC=NO");
144 : }
145 16 : action->readInputLine( labout + "_lessthan: SUM ARG=" + sum_arg + " PERIODIC=NO");
146 : }
147 294 : if( keymap.count("LESS_THAN1") ) {
148 2 : for(unsigned i=1;; ++i) {
149 : std::string istr;
150 6 : Tools::convert( i, istr );
151 12 : if( !keymap.count("LESS_THAN" + istr ) ) {
152 : break;
153 : }
154 12 : std::string sum_arg = labout + "_lt" + istr, lt_string1 = keymap.find("LESS_THAN" + istr)->second;
155 8 : action->readInputLine( labout + "_lt" + istr + ": LESS_THAN ARG=" + argin + " SWITCH={" + lt_string1 + "}");
156 4 : if( weights.length()>0 ) {
157 0 : sum_arg = labout + "_wlt" + istr;
158 0 : action->readInputLine( labout + "_wlt" + istr + ": CUSTOM ARG=" + weights + "," + labout + "_lt" + istr + " FUNC=x*y PERIODIC=NO");
159 : }
160 8 : action->readInputLine( labout + "_lessthan-" + istr + ": SUM ARG=" + sum_arg + " PERIODIC=NO");
161 4 : }
162 : }
163 : // Parse MORE_THAN
164 294 : if( keymap.count("MORE_THAN") ) {
165 48 : std::string sum_arg=labout + "_mt", mt_string = keymap.find("MORE_THAN")->second;
166 48 : action->readInputLine( labout + "_mt: MORE_THAN ARG=" + argin + " SWITCH={" + mt_string + "}");
167 24 : if( weights.length()>0 ) {
168 1 : sum_arg = labout + "_wmt";
169 2 : action->readInputLine( labout + "_wmt: CUSTOM ARG=" + weights + "," + labout + "_mt FUNC=x*y PERIODIC=NO" );
170 : }
171 48 : action->readInputLine( labout + "_morethan: SUM ARG=" + sum_arg + " PERIODIC=NO");
172 : }
173 294 : if( keymap.count("MORE_THAN1") ) {
174 0 : for(unsigned i=1;; ++i) {
175 : std::string istr;
176 0 : Tools::convert( i, istr );
177 0 : if( !keymap.count("MORE_THAN" + istr ) ) {
178 : break;
179 : }
180 0 : std::string sum_arg = labout + "_mt" + istr, mt_string1 = keymap.find("MORE_THAN" + istr)->second;
181 0 : action->readInputLine( labout + "_mt" + istr + ": MORE_THAN ARG=" + argin + " SWITCH={" + mt_string1 + "}");
182 0 : if( weights.length()>0 ) {
183 0 : sum_arg = labout + "_wmt" + istr;
184 0 : action->readInputLine( labout + "_wmt" + istr + ": CUSTOM ARG=" + weights + "," + labout + "_lt" + istr + " FUNC=x*y PERIODIC=NO");
185 : }
186 0 : action->readInputLine( labout + "_morethan-" + istr + ": SUM ARG=" + sum_arg + " PERIODIC=NO");
187 0 : }
188 : }
189 : // Parse ALT_MIN
190 294 : if( keymap.count("ALT_MIN") ) {
191 1 : if( weights.length()>0 ) {
192 0 : plumed_merror("cannot use ALT_MIN with this shortcut");
193 : }
194 2 : std::string amin_string = keymap.find("ALT_MIN")->second;
195 1 : std::size_t dd = amin_string.find("BETA");
196 1 : std::string beta_str = amin_string.substr(dd+5);
197 1 : beta_str.erase(std::remove_if(beta_str.begin(), beta_str.end(), ::isspace), beta_str.end());
198 2 : action->readInputLine( labout + "_me_altmin: CUSTOM ARG=" + argin + " FUNC=exp(-x*" + beta_str + ") PERIODIC=NO");
199 2 : action->readInputLine( labout + "_mec_altmin: SUM ARG=" + labout + "_me_altmin PERIODIC=NO");
200 2 : action->readInputLine( labout + "_altmin: CUSTOM ARG=" + labout + "_mec_altmin FUNC=-log(x)/" + beta_str + " PERIODIC=NO");
201 : }
202 : // Parse MIN
203 294 : if( keymap.count("MIN") ) {
204 1 : if( weights.length()>0 ) {
205 0 : plumed_merror("cannot use MIN with this shortcut");
206 : }
207 2 : std::string min_string = keymap.find("MIN")->second;
208 1 : std::size_t dd = min_string.find("BETA");
209 1 : std::string beta_str = min_string.substr(dd+5);
210 1 : beta_str.erase(std::remove_if(beta_str.begin(), beta_str.end(), ::isspace), beta_str.end());
211 2 : action->readInputLine( labout + "_me_min: CUSTOM ARG=" + argin + " FUNC=exp(" + beta_str + "/x) PERIODIC=NO");
212 2 : action->readInputLine( labout + "_mec_min: SUM ARG=" + labout + "_me_min PERIODIC=NO");
213 2 : action->readInputLine( labout + "_min: CUSTOM ARG=" + labout + "_mec_min FUNC=" + beta_str + "/log(x) PERIODIC=NO");
214 : }
215 : // Parse MAX
216 294 : if( keymap.count("MAX") ) {
217 3 : if( weights.length()>0 ) {
218 0 : plumed_merror("cannot use MAX with this shortcut");
219 : }
220 6 : std::string max_string = keymap.find("MAX")->second;
221 3 : std::size_t dd = max_string.find("BETA");
222 3 : std::string beta_str = max_string.substr(dd+5);
223 3 : beta_str.erase(std::remove_if(beta_str.begin(), beta_str.end(), ::isspace), beta_str.end());
224 6 : action->readInputLine( labout + "_me_max: CUSTOM ARG=" + argin + " FUNC=exp(x/" + beta_str + ") PERIODIC=NO");
225 6 : action->readInputLine( labout + "_mec_max: SUM ARG=" + labout + "_me_max PERIODIC=NO");
226 6 : action->readInputLine( labout + "_max: CUSTOM ARG=" + labout + "_mec_max FUNC=" + beta_str + "*log(x) PERIODIC=NO");
227 : }
228 : // Parse HIGHEST
229 294 : if( keymap.count("HIGHEST") ) {
230 8 : if( weights.length()>0 ) {
231 0 : plumed_merror("cannot use HIGHEST with this shortcut");
232 : }
233 16 : action->readInputLine( labout + "_highest: HIGHEST ARG=" + argin );
234 : }
235 : // Parse LOWEST
236 294 : if( keymap.count("LOWEST") ) {
237 3 : if( weights.length()>0 ) {
238 0 : plumed_merror("cannot use LOWEST with this shortcut");
239 : }
240 6 : action->readInputLine( labout + "_lowest: LOWEST ARG=" + argin );
241 : }
242 : // Parse SUM
243 294 : if( keymap.count("SUM") ) {
244 44 : std::string sum_arg=argin;
245 44 : if( weights.length()>0 ) {
246 22 : sum_arg = labout + "_wsum";
247 44 : action->readInputLine( labout + "_wsum: CUSTOM ARG=" + weights + "," + argin + " FUNC=x*y PERIODIC=NO");
248 : }
249 88 : action->readInputLine( labout + "_sum: SUM ARG=" + sum_arg + " PERIODIC=NO");
250 : }
251 : // Parse MEAN
252 294 : if( keymap.count("MEAN") ) {
253 64 : if( weights.length()>0 ) {
254 0 : plumed_merror("cannot use MEAN with this shortcut");
255 : }
256 128 : action->readInputLine( labout + "_mean: MEAN ARG=" + argin + " PERIODIC=NO");
257 : }
258 : // Parse BETWEEN
259 294 : if( keymap.count("BETWEEN") ) {
260 8 : std::string sum_arg=labout + "_bt";
261 16 : std::string bt_string = keymap.find("BETWEEN")->second;
262 16 : action->readInputLine( labout + "_bt: BETWEEN ARG=" + argin + " SWITCH={" + bt_string + "}" );
263 8 : if( weights.length()>0 ) {
264 0 : sum_arg = labout + "_wbt";
265 0 : action->readInputLine( labout + "_wbt: CUSTOM ARG=" + weights + "," + labout + "_bt FUNC=x*y PERIODIC=NO");
266 : }
267 16 : action->readInputLine( labout + "_between: SUM ARG=" + sum_arg + " PERIODIC=NO");
268 : }
269 :
270 294 : if( keymap.count("BETWEEN1") ) {
271 6 : for(unsigned i=1;; ++i) {
272 : std::string istr;
273 17 : Tools::convert( i, istr );
274 34 : if( !keymap.count("BETWEEN" + istr) ) {
275 : break;
276 : }
277 11 : std::string sum_arg=labout + "_bt" + istr;
278 22 : std::string bt_string1 = keymap.find("BETWEEN" + istr)->second;
279 22 : action->readInputLine( labout + "_bt" + istr + ": BETWEEN ARG=" + argin + " SWITCH={" + bt_string1 + "}" );
280 11 : if( weights.length()>0 ) {
281 2 : sum_arg = labout + "_wbt" + istr;
282 2 : action->readInputLine( labout + "_wbt" + istr + ": CUSTOM ARG=" + weights + "," + labout + "_bt" + istr + " FUNC=x*y PERIODIC=NO");
283 : }
284 22 : action->readInputLine( labout + "_between-" + istr + ": SUM ARG=" + sum_arg + " PERIODIC=NO");
285 11 : }
286 : }
287 : // Parse HISTOGRAM
288 294 : if( keymap.count("HISTOGRAM") ) {
289 24 : std::vector<std::string> words=Tools::getWords( keymap.find("HISTOGRAM")->second );
290 : unsigned nbins;
291 12 : bool found=Tools::parse(words,"NBINS",nbins,0); // Need replica index
292 12 : if( !found ) {
293 0 : plumed_merror("did not find NBINS in specification for HISTOGRAM");
294 : }
295 : double lower;
296 12 : found=Tools::parse(words,"LOWER",lower,0);
297 12 : if( !found ) {
298 0 : plumed_merror("did not find LOWER in specification for HISTOGRAM");
299 : }
300 : double upper;
301 12 : found=Tools::parse(words,"UPPER",upper,0);
302 12 : if( !found ) {
303 0 : plumed_merror("did not find UPPER in specification for HISTOGRAM");
304 : }
305 12 : double delr = ( upper - lower ) / static_cast<double>( nbins );
306 12 : double smear=0.5;
307 12 : found=Tools::parse(words,"SMEAR",smear,0);
308 12 : if( !found ) {
309 5 : smear = 0.5;
310 : }
311 41 : for(unsigned i=0; i<nbins; ++i) {
312 : std::string smstr, istr;
313 29 : Tools::convert( i+1, istr );
314 29 : Tools::convert( smear, smstr );
315 58 : std::string sum_arg=labout + "_bt" + istr;
316 : std::string low_str, high_str;
317 29 : Tools::convert( lower + i*delr, low_str );
318 29 : Tools::convert( lower + (i+1)*delr, high_str );
319 58 : action->readInputLine( labout + "_bt" + istr + ": BETWEEN ARG=" + argin + " SWITCH={" + words[0] + " LOWER=" + low_str + " UPPER=" + high_str + " SMEAR=" + smstr + "}");
320 29 : if( weights.length()>0 ) {
321 0 : sum_arg = labout + "_wbt" + istr;
322 0 : action->readInputLine( labout + "_wbt" + istr + ": CUSTOM ARG=" + weights + "," + labout + "_bt" + istr + " FUNC=x*y PERIODIC=NO");
323 : }
324 58 : action->readInputLine( labout + "_between-" + istr + ": SUM ARG=" + sum_arg + " PERIODIC=NO");
325 : }
326 12 : }
327 : }
328 :
329 : }
330 : }
|