Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-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 "ActionWithArguments.h"
23 : #include "ActionWithValue.h"
24 : #include "tools/PDB.h"
25 : #include "PlumedMain.h"
26 : #include "ActionSet.h"
27 : #include <iostream>
28 : #ifdef __PLUMED_HAS_CREGEX
29 : #include <cstring>
30 : #include <regex.h>
31 : #endif
32 :
33 : namespace PLMD {
34 :
35 3409 : void ActionWithArguments::registerKeywords(Keywords& keys) {
36 6818 : keys.reserve("numbered","ARG","the input for this action is the scalar output from one or more other actions. The particular scalars that you will use "
37 : "are referenced using the label of the action. If the label appears on its own then it is assumed that the Action calculates "
38 : "a single scalar value. The value of this scalar is thus used as the input to this new action. If * or *.* appears the "
39 : "scalars calculated by all the proceeding actions in the input file are taken. Some actions have multi-component outputs and "
40 : "each component of the output has a specific label. For example a \\ref DISTANCE action labelled dist may have three components "
41 : "x, y and z. To take just the x component you should use dist.x, if you wish to take all three components then use dist.*."
42 : "More information on the referencing of Actions can be found in the section of the manual on the PLUMED \\ref Syntax. "
43 : "Scalar values can also be "
44 : "referenced using POSIX regular expressions as detailed in the section on \\ref Regex. To use this feature you you must compile "
45 : "PLUMED with the appropriate flag.");
46 3409 : }
47 :
48 3065 : void ActionWithArguments::parseArgumentList(const std::string&key,std::vector<Value*>&arg) {
49 : std::string def;
50 : std::vector<std::string> c;
51 : arg.clear();
52 3065 : parseVector(key,c);
53 3645 : if( c.size()==0 && (keywords.style(key,"compulsory") || keywords.style(key,"hidden")) ) {
54 0 : if( keywords.getDefaultValue(key,def) ) {
55 0 : c.push_back( def );
56 : } else {
57 : return;
58 : }
59 : }
60 3065 : interpretArgumentList(c,arg);
61 3065 : }
62 :
63 56 : bool ActionWithArguments::parseArgumentList(const std::string&key,int i,std::vector<Value*>&arg) {
64 : std::vector<std::string> c;
65 : arg.clear();
66 56 : if(parseNumberedVector(key,i,c)) {
67 14 : interpretArgumentList(c,arg);
68 : return true;
69 : } else {
70 : return false;
71 : }
72 56 : }
73 :
74 3175 : void ActionWithArguments::interpretArgumentList(const std::vector<std::string>& c, std::vector<Value*>&arg) {
75 8886 : for(unsigned i=0; i<c.size(); i++) {
76 : // is a regex? then just interpret it. The signal is ()
77 5714 : if(!c[i].compare(0,1,"(")) {
78 201 : unsigned l=c[i].length();
79 201 : if(!c[i].compare(l-1,1,")")) {
80 : // start regex parsing
81 : #ifdef __PLUMED_HAS_CREGEX
82 : // take the string enclosed in quotes and put in round brackets
83 200 : std::string myregex=c[i];
84 : //log<<" Evaluating regexp for this action: "<<myregex<<"\n";
85 :
86 : regex_t reg; // regular expression
87 :
88 200 : int errcode=regcomp(®, myregex.c_str(),REG_EXTENDED|REG_NEWLINE); // compile the regular expression
89 200 : if(errcode) {
90 : // one can check the errors asking to regerror
91 1 : size_t errbuf_size = regerror(errcode, ®, NULL, 0);
92 2 : std::vector<char> errbuf(errbuf_size);
93 1 : regerror(errcode, ®, errbuf.data(), errbuf_size);
94 3 : plumed_error()<<"Error parsing regular expression: "<<errbuf.data();
95 : }
96 :
97 : // call regfree when reg goes out of scope
98 : auto deleter=[](regex_t* r) {
99 199 : regfree(r);
100 : };
101 : std::unique_ptr<regex_t,decltype(deleter)> reg_deleter(®,deleter);
102 :
103 199 : plumed_massert(reg.re_nsub==1,"I can parse with only one subexpression");
104 : regmatch_t match;
105 : // select all the actions that have a value
106 199 : std::vector<ActionWithValue*> all=plumed.getActionSet().select<ActionWithValue*>();
107 199 : if( all.empty() ) {
108 0 : error("your input file is not telling plumed to calculate anything");
109 : }
110 : bool found_something=false;
111 1522 : for(unsigned j=0; j<all.size(); j++) {
112 1323 : std::vector<std::string> ss=all[j]->getComponentsVector();
113 220386 : for(unsigned k=0; k<ss.size(); ++k) {
114 219063 : unsigned ll=std::strlen(ss[k].c_str())+1;
115 219063 : std::vector<char> str(ll);
116 : std::strcpy(&str[0],ss[k].c_str());
117 : const char *ppstr=&str[0];
118 219063 : if(!regexec(®, ppstr, reg.re_nsub, &match, 0)) {
119 : //log.printf(" Something matched with \"%s\" : ",ss[k].c_str());
120 : do {
121 27856 : if (match.rm_so != -1) { /* The regex is matching part of a string */
122 27856 : size_t matchlen = match.rm_eo - match.rm_so;
123 27856 : std::vector<char> submatch(matchlen+1);
124 27856 : std::strncpy(submatch.data(), ppstr+match.rm_so, matchlen+1);
125 27856 : submatch[matchlen]='\0';
126 : //log.printf(" subpattern %s\n", submatch.data());
127 : // this is the match: try to see if it is a valid action
128 27856 : std::string putativeVal(submatch.data());
129 27856 : if( all[j]->exists(putativeVal) ) {
130 21369 : arg.push_back(all[j]->copyOutput(putativeVal));
131 : found_something=true;
132 : //log.printf(" Action %s added! \n",putativeVal.c_str());
133 : }
134 : }
135 27856 : ppstr += match.rm_eo; /* Restart from last match */
136 27856 : } while(!regexec(®,ppstr,reg.re_nsub,&match,0));
137 : }
138 : }
139 1323 : }
140 199 : if(!found_something) {
141 0 : plumed_error()<<"There isn't any action matching your regex " << myregex;
142 : }
143 : #else
144 : plumed_merror("Regexp support not compiled!");
145 : #endif
146 : } else {
147 2 : plumed_merror("did you want to use regexp to input arguments? enclose it between two round braces (...) with no spaces!");
148 : }
149 : } else {
150 : std::size_t dot=c[i].find_first_of('.');
151 5513 : std::string a=c[i].substr(0,dot);
152 5513 : std::string name=c[i].substr(dot+1);
153 5513 : if(c[i].find(".")!=std::string::npos) { // if it contains a dot:
154 2046 : if(a=="*" && name=="*") {
155 : // Take all values from all actions
156 1 : std::vector<ActionWithValue*> all=plumed.getActionSet().select<ActionWithValue*>();
157 1 : if( all.empty() ) {
158 0 : error("your input file is not telling plumed to calculate anything");
159 : }
160 9 : for(unsigned j=0; j<all.size(); j++) {
161 18 : for(int k=0; k<all[j]->getNumberOfComponents(); ++k) {
162 10 : arg.push_back(all[j]->copyOutput(k));
163 : }
164 : }
165 2029 : } else if ( name=="*") {
166 : // Take all the values from an action with a specific name
167 586 : ActionWithValue* action=plumed.getActionSet().selectWithLabel<ActionWithValue*>(a);
168 586 : if(!action) {
169 0 : std::string str=" (hint! the actions with value in this ActionSet are: ";
170 0 : str+=plumed.getActionSet().getLabelList<ActionWithValue*>()+")";
171 0 : error("cannot find action named " + a + str);
172 : }
173 586 : if( action->getNumberOfComponents()==0 ) {
174 0 : error("found " + a +".* indicating use all components calculated by action with label " + a + " but this action has no components");
175 : }
176 5791 : for(int k=0; k<action->getNumberOfComponents(); ++k) {
177 5205 : arg.push_back(action->copyOutput(k));
178 : }
179 1443 : } else if ( a=="*" ) {
180 : // Take components from all actions with a specific name
181 15 : std::vector<ActionWithValue*> all=plumed.getActionSet().select<ActionWithValue*>();
182 15 : if( all.empty() ) {
183 0 : error("your input file is not telling plumed to calculate anything");
184 : }
185 : unsigned nval=0;
186 78 : for(unsigned j=0; j<all.size(); j++) {
187 : std::string flab;
188 126 : flab=all[j]->getLabel() + "." + name;
189 63 : if( all[j]->exists(flab) ) {
190 42 : arg.push_back(all[j]->copyOutput(flab));
191 42 : nval++;
192 : }
193 : }
194 15 : if(nval==0) {
195 0 : error("found no actions with a component called " + name );
196 : }
197 : } else {
198 : // Take values with a specific name
199 1428 : ActionWithValue* action=plumed.getActionSet().selectWithLabel<ActionWithValue*>(a);
200 1428 : if(!action) {
201 0 : std::string str=" (hint! the actions with value in this ActionSet are: ";
202 0 : str+=plumed.getActionSet().getLabelList<ActionWithValue*>()+")";
203 0 : error("cannot find action named " + a +str);
204 : }
205 1428 : if( !(action->exists(c[i])) ) {
206 0 : std::string str=" (hint! the components in this actions are: ";
207 0 : str+=action->getComponentsList()+")";
208 0 : error("action " + a + " has no component named " + name + str);
209 : } ;
210 1428 : arg.push_back(action->copyOutput(c[i]));
211 : }
212 : } else { // if it doesn't contain a dot
213 3483 : if(c[i]=="*") {
214 : // Take all values from all actions
215 108 : std::vector<ActionWithValue*> all=plumed.getActionSet().select<ActionWithValue*>();
216 108 : if( all.empty() ) {
217 0 : error("your input file is not telling plumed to calculate anything");
218 : }
219 732 : for(unsigned j=0; j<all.size(); j++) {
220 1340 : for(int k=0; k<all[j]->getNumberOfComponents(); ++k) {
221 716 : arg.push_back(all[j]->copyOutput(k));
222 : }
223 : }
224 : } else {
225 3375 : ActionWithValue* action=plumed.getActionSet().selectWithLabel<ActionWithValue*>(c[i]);
226 3375 : if(!action) {
227 1 : std::string str=" (hint! the actions with value in this ActionSet are: ";
228 2 : str+=plumed.getActionSet().getLabelList<ActionWithValue*>()+")";
229 3 : error("cannot find action named " + c[i] + str );
230 : }
231 3374 : if( !(action->exists(c[i])) ) {
232 0 : std::string str=" (hint! the components in this actions are: ";
233 0 : str+=action->getComponentsList()+")";
234 0 : error("action " + c[i] + " has no component named " + c[i] +str);
235 : };
236 3375 : arg.push_back(action->copyOutput(c[i]));
237 : }
238 : }
239 : }
240 : }
241 3172 : }
242 :
243 416 : void ActionWithArguments::expandArgKeywordInPDB( const PDB& pdb ) {
244 416 : std::vector<std::string> arg_names = pdb.getArgumentNames();
245 416 : if( arg_names.size()>0 ) {
246 : std::vector<Value*> arg_vals;
247 80 : interpretArgumentList( arg_names, arg_vals );
248 : }
249 416 : }
250 :
251 179128 : void ActionWithArguments::requestArguments(const std::vector<Value*> &arg) {
252 179128 : plumed_massert(!lockRequestArguments,"requested argument list can only be changed in the prepare() method");
253 179128 : arguments=arg;
254 179128 : clearDependencies();
255 : std::string fullname;
256 : std::string name;
257 1255830 : for(unsigned i=0; i<arguments.size(); i++) {
258 1076702 : fullname=arguments[i]->getName();
259 1076702 : if(fullname.find(".")!=std::string::npos) {
260 : std::size_t dot=fullname.find_first_of('.');
261 739312 : name=fullname.substr(0,dot);
262 : } else {
263 : name=fullname;
264 : }
265 1076702 : ActionWithValue* action=plumed.getActionSet().selectWithLabel<ActionWithValue*>(name);
266 1076702 : plumed_massert(action,"cannot find action named (in requestArguments - this is weird)" + name);
267 1076702 : addDependency(action);
268 : }
269 179128 : }
270 :
271 4 : void ActionWithArguments::requestExtraDependencies(const std::vector<Value*> &extra) {
272 4 : plumed_massert(!lockRequestArguments,"requested argument list can only be changed in the prepare() method");
273 : std::string fullname;
274 : std::string name;
275 9 : for(unsigned i=0; i<extra.size(); i++) {
276 5 : fullname=extra[i]->getName();
277 5 : if(fullname.find(".")!=std::string::npos) {
278 : std::size_t dot=fullname.find_first_of('.');
279 0 : name=fullname.substr(0,dot);
280 : } else {
281 : name=fullname;
282 : }
283 5 : ActionWithValue* action=plumed.getActionSet().selectWithLabel<ActionWithValue*>(name);
284 5 : plumed_massert(action,"cannot find action named (in requestArguments - this is weird)" + name);
285 5 : addDependency(action);
286 : }
287 4 : }
288 :
289 3012 : ActionWithArguments::ActionWithArguments(const ActionOptions&ao):
290 : Action(ao),
291 3012 : lockRequestArguments(false) {
292 6027 : if( keywords.exists("ARG") ) {
293 : std::vector<Value*> arg;
294 5838 : parseArgumentList("ARG",arg);
295 :
296 2916 : if(!arg.empty()) {
297 2739 : log.printf(" with arguments");
298 26421 : for(unsigned i=0; i<arg.size(); i++) {
299 23682 : log.printf(" %s",arg[i]->getName().c_str());
300 : }
301 2739 : log.printf("\n");
302 : }
303 2916 : requestArguments(arg);
304 : }
305 3009 : }
306 :
307 58 : void ActionWithArguments::calculateNumericalDerivatives( ActionWithValue* a ) {
308 58 : if(!a) {
309 58 : a=dynamic_cast<ActionWithValue*>(this);
310 58 : plumed_massert(a,"cannot compute numerical derivatives for an action without values");
311 : }
312 :
313 58 : const size_t nval=a->getNumberOfComponents();
314 : const size_t npar=arguments.size();
315 58 : std::vector<double> value (nval*npar);
316 161 : for(int i=0; i<npar; i++) {
317 103 : double arg0=arguments[i]->get();
318 103 : arguments[i]->set(arg0+std::sqrt(epsilon));
319 103 : a->calculate();
320 103 : arguments[i]->set(arg0);
321 1367 : for(int j=0; j<nval; j++) {
322 1264 : value[i*nval+j]=a->getOutputQuantity(j);
323 : }
324 : }
325 58 : a->calculate();
326 58 : a->clearDerivatives();
327 1192 : for(int j=0; j<nval; j++) {
328 1134 : Value* v=a->copyOutput(j);
329 1134 : if( v->hasDerivatives() )
330 670 : for(int i=0; i<npar; i++) {
331 400 : v->addDerivative(i,(value[i*nval+j]-a->getOutputQuantity(j))/std::sqrt(epsilon));
332 : }
333 : }
334 58 : }
335 :
336 261 : double ActionWithArguments::getProjection(unsigned i,unsigned j)const {
337 261 : plumed_massert(i<arguments.size()," making projections with an index which is too large");
338 261 : plumed_massert(j<arguments.size()," making projections with an index which is too large");
339 261 : const Value* v1=arguments[i];
340 261 : const Value* v2=arguments[j];
341 261 : return Value::projection(*v1,*v2);
342 : }
343 :
344 482 : void ActionWithArguments::addForcesOnArguments( const std::vector<double>& forces ) {
345 606 : for(unsigned i=0; i<arguments.size(); ++i) {
346 124 : arguments[i]->addForce( forces[i] );
347 : }
348 482 : }
349 :
350 : }
|
