Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2012-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 "CLTool.h"
23 : #include "core/CLToolRegister.h"
24 : #include "tools/Tools.h"
25 : #include "core/Action.h"
26 : #include "core/ActionRegister.h"
27 : #include "core/PlumedMain.h"
28 : #include "tools/Communicator.h"
29 : #include "tools/Random.h"
30 : #include <cstdio>
31 : #include <string>
32 : #include <vector>
33 : #include "tools/File.h"
34 : #include "core/Value.h"
35 : #include "tools/Matrix.h"
36 :
37 : namespace PLMD {
38 : namespace cltools {
39 :
40 : //+PLUMEDOC TOOLS sum_hills
41 : /*
42 : sum_hills is a tool that allows one to to use plumed to post-process an existing hills/colvar file
43 :
44 : This tool is most frequently used as shown below after a [METAD](METAD.md) simulation has been performed to sum all the Gaussians in the hills file
45 : and output the free energy surface.
46 :
47 : ```plumed
48 : plumed sum_hills --hills PATHTOMYHILLSFILE
49 : ```
50 :
51 : The default name for the output file will be `fes.dat`
52 : Note that plumed automatically detects the
53 : number of the variables you have and their periodicity.
54 : Additionally, if you use flexible hills (multivariate Gaussian kernels), plumed will understand it from the HILLS file.
55 :
56 : The sum_hills tool can also accept multiple files that will be integrated one after the other as shown below
57 :
58 : ```plumed
59 : plumed sum_hills --hills PATHTOMYHILLSFILE1,PATHTOMYHILLSFILE2,PATHTOMYHILLSFILE3
60 : ```
61 :
62 : if you want to integrate out some variable from your output free energy surface you do
63 :
64 : ```plumed
65 : plumed sum_hills --hills PATHTOMYHILLSFILE --idw t1 --kt 0.6
66 : ```
67 :
68 : where with `--idw` you define the variables that you want in the output free energy surface.
69 : All other variables in the hills file will be integrated out. `--kt` defines the temperature of the system in energy units.
70 : (be consistent with the units you have in your hills: plumed will not check this for you)
71 : If you need more variables then you may use a comma separated syntax shown below:
72 :
73 : ```plumed
74 : plumed sum_hills --hills PATHTOMYHILLSFILE --idw t1,t2 --kt 0.6
75 : ```
76 :
77 : You can define the output grid with the number of bins you want by using a command like this one
78 :
79 : ```plumed
80 : plumed sum_hills --bin 99,99 --hills PATHTOMYHILLSFILE
81 : ```
82 :
83 : In the command above the min/max to use are detected for you. If you want to specify the min/max yourself you can do so by using
84 : a command like the one shown below:
85 :
86 : ```plumed
87 : plumed sum_hills --bin 99,99 --min -pi,-pi --max pi,pi --hills PATHTOMYHILLSFILE
88 : ```
89 :
90 : You can of course use numbers instead of -pi/pi.
91 :
92 : You can use a `--stride` keyword to have a dump of the free energy estimate after each bunch of hills you read as shown in the following command:
93 :
94 : ```plumed
95 : plumed sum_hills --stride 300 --hills PATHTOMYHILLSFILE
96 : ```
97 :
98 : If you have run well tempered metadynamics you can also ask `sum_hills` to output the negative bias instead of the free energy by using the keyword `--negbias`
99 : as shown below
100 :
101 : ```plumed
102 : plumed sum_hills --negbias --hills PATHTOMYHILLSFILE
103 : ```
104 :
105 : Here the default output file name will be negativebias.dat
106 :
107 : From time to time you might need to use HILLS or a COLVAR file
108 : as it was just a simple set of points from which you want to build
109 : a free energy by using -(1/beta)log(P). If you want to do this operation then
110 : you use the `--histo` option as shown below:
111 :
112 : ```plumed
113 : plumed sum_hills --histo PATHTOMYCOLVARORHILLSFILE --sigma 0.2,0.2 --kt 0.6
114 : ```
115 :
116 : in this case you need a `--kt` to do the reweighting and you
117 : need to set bandwidth parameters using the `--sigma` option for the histogram calculation as this histogram is computed using
118 : Gaussian kernels, so it will be a continuous histogram.
119 :
120 : For the command above the default output is called histo.dat.
121 : Note that also here you can have multiple input files separated by a comma.
122 :
123 : Additionally, if you want to compute the histogram and sum of the hills from the same file you can do this by using a command like this one:
124 :
125 : ```plumed
126 : plumed sum_hills --hills --histo PATHTOMYCOLVARORHILLSFILE --sigma 0.2,0.2 --kt 0.6
127 : ```
128 :
129 : The two files can be eventually the same
130 :
131 : Another interesting thing one can do is monitor the difference in blocks as a metadynamics goes on.
132 : When the bias deposited is constant over the whole domain one can consider to be at convergence.
133 : This can be done with the `--nohistory` keyword
134 :
135 : ```plumed
136 : plumed sum_hills --stride 300 --hills PATHTOMYHILLSFILE --nohistory
137 : ```
138 :
139 : and similarly one can do the same for an histogram file
140 :
141 : ```plumed
142 : plumed sum_hills --histo PATHTOMYCOLVARORHILLSFILE --sigma 0.2,0.2 --kt 0.6 --nohistory
143 : ```
144 :
145 : just to check the hypothetical free energy calculated in single blocks of time during a simulation
146 : and not in a cumulative way
147 :
148 : The output format can be controlled by using the `--fmt` option as shown below
149 :
150 : ```plumed
151 : plumed sum_hills --hills PATHTOMYHILLSFILE --fmt %8.3f
152 : ```
153 :
154 : where here we chose a float with length of 8 and 3 digits
155 :
156 : The output can be named in a arbitrary way by using the `--output` option as shown below:
157 :
158 : ```plumed
159 : plumed sum_hills --hills PATHTOMYHILLSFILE --outfile myfes.dat
160 : ```
161 :
162 : This command produces a file myfes.dat which contains the free energy surface.
163 :
164 : If you use stride, the name specied with `--output` is the suffix so this command:
165 :
166 : ```plumed
167 : plumed sum_hills --hills PATHTOMYHILLSFILE --outfile myfes_ --stride 100
168 : ```
169 :
170 : produces output files called myfes_0.dat, myfes_1.dat, myfes_2.dat etc.
171 :
172 : The same is true for the output coming from histogram that is used. So this example
173 :
174 : ```plumed
175 : plumed sum_hills --histo HILLS --kt 2.5 --sigma 0.01 --outhisto myhisto.dat
176 : ```
177 :
178 : produces a file myhisto.dat, while this input
179 :
180 : ```plumed
181 : plumed sum_hills --histo HILLS --kt 2.5 --sigma 0.01 --outhisto myhisto_ --stride 100
182 : ```
183 :
184 : produces output files called `myhisto_0.dat`, `myhisto_1.dat`, `myhisto_3.dat` etc..
185 :
186 : */
187 : //+ENDPLUMEDOC
188 :
189 : class CLToolSumHills : public CLTool {
190 : public:
191 : static void registerKeywords( Keywords& keys );
192 : explicit CLToolSumHills(const CLToolOptions& co );
193 : int main(FILE* in,FILE*out,Communicator& pc) override;
194 : std::string description()const override;
195 : /// find a list of variables present, if they are periodic and which is the period
196 : /// return false if the file does not exist
197 : static bool findCvsAndPeriodic(const std::string & filename, std::vector< std::vector <std::string> > &cvs,std::vector<std::string> &pmin,std::vector<std::string> &pmax, bool &multivariate, std::string &lowI_, std::string &uppI_);
198 : };
199 :
200 5683 : void CLToolSumHills::registerKeywords( Keywords& keys ) {
201 5683 : CLTool::registerKeywords( keys );
202 5683 : keys.addFlag("--help-debug",false,"print special options that can be used to create regtests");
203 5683 : keys.add("optional","--hills","specify the name of the hills file");
204 5683 : keys.add("optional","--histo","specify the name of the file for histogram a colvar/hills file is good");
205 5683 : keys.add("optional","--stride","specify the stride for integrating hills file (default 0=never)");
206 5683 : keys.add("optional","--min","the lower bounds for the grid");
207 5683 : keys.add("optional","--max","the upper bounds for the grid");
208 5683 : keys.add("optional","--bin","the number of bins for the grid");
209 5683 : keys.add("optional","--spacing","grid spacing, alternative to the number of bins");
210 5683 : keys.add("optional","--idw","specify the variables to be used for the free-energy/histogram (default is all). With --hills the other variables will be integrated out, with --histo the other variables won't be considered");
211 5683 : keys.add("optional","--outfile","specify the output file for sumhills");
212 5683 : keys.add("optional","--outhisto","specify the output file for the histogram");
213 5683 : keys.add("optional","--kt","specify temperature in energy units for integrating out variables");
214 5683 : keys.add("optional","--sigma"," a vector that specify the sigma for binning (only needed when doing histogram ");
215 5683 : keys.addFlag("--negbias",false," print the negative bias instead of the free energy (only needed with well tempered runs and flexible hills) ");
216 5683 : keys.addFlag("--nohistory",false," to be used with --stride: it splits the bias/histogram in pieces without previous history ");
217 5683 : keys.addFlag("--mintozero",false," it translate all the minimum value in bias/histogram to zero (useful to compare results) ");
218 5683 : keys.add("optional","--fmt","specify the output format");
219 5683 : }
220 :
221 14 : CLToolSumHills::CLToolSumHills(const CLToolOptions& co ):
222 14 : CLTool(co) {
223 14 : inputdata=inputType::commandline;
224 14 : }
225 :
226 5 : std::string CLToolSumHills::description()const {
227 5 : return "sum the hills with plumed";
228 : }
229 :
230 9 : int CLToolSumHills::main(FILE* in,FILE*out,Communicator& pc) {
231 :
232 : // Read the hills input file name
233 : std::vector<std::string> hillsFiles;
234 : bool dohills;
235 18 : dohills=parseVector("--hills",hillsFiles);
236 : // Read the histogram file
237 : std::vector<std::string> histoFiles;
238 : bool dohisto;
239 9 : dohisto=parseVector("--histo",histoFiles);
240 :
241 9 : plumed_massert(dohisto || dohills,"you should use --histo or/and --hills command");
242 :
243 : std::vector< std::vector<std::string> > vcvs;
244 : std::vector<std::string> vpmin;
245 : std::vector<std::string> vpmax;
246 : std::string lowI_, uppI_;
247 9 : if(dohills) {
248 : // parse it as it was a restart
249 : bool vmultivariate;
250 8 : findCvsAndPeriodic(hillsFiles[0], vcvs, vpmin, vpmax, vmultivariate, lowI_, uppI_);
251 : }
252 :
253 : std::vector< std::vector<std::string> > hcvs;
254 : std::vector<std::string> hpmin;
255 : std::vector<std::string> hpmax;
256 :
257 : std::vector<std::string> sigma;
258 9 : if(dohisto) {
259 : bool hmultivariate;
260 1 : findCvsAndPeriodic(histoFiles[0], hcvs, hpmin, hpmax, hmultivariate, lowI_, uppI_);
261 : // here need also the vector of sigmas
262 2 : parseVector("--sigma",sigma);
263 1 : if(sigma.size()==0) {
264 0 : plumed_merror("you should define --sigma vector when using histogram");
265 : }
266 : lowI_=uppI_="-1."; // Interval is not use for histograms
267 : }
268 :
269 9 : if(dohisto && dohills) {
270 0 : plumed_massert(vcvs==hcvs,"variables for histogram and bias should have the same labels");
271 0 : plumed_massert(hpmin==vpmin,"variables for histogram and bias should have the same min for periodicity");
272 0 : plumed_massert(hpmax==vpmax,"variables for histogram and bias should have the same max for periodicity");
273 : }
274 :
275 : // now put into a neutral vector
276 :
277 : std::vector< std::vector<std::string> > cvs;
278 : std::vector<std::string> pmin;
279 : std::vector<std::string> pmax;
280 :
281 9 : if(dohills) {
282 8 : cvs=vcvs;
283 8 : pmin=vpmin;
284 8 : pmax=vpmax;
285 : }
286 9 : if(dohisto) {
287 1 : cvs=hcvs;
288 1 : pmin=hpmin;
289 1 : pmax=hpmax;
290 : }
291 :
292 :
293 : // setup grids
294 : unsigned grid_check=0;
295 9 : std::vector<std::string> gmin(cvs.size());
296 18 : if(parseVector("--min",gmin)) {
297 4 : if(gmin.size()!=cvs.size() && gmin.size()!=0) {
298 0 : plumed_merror("not enough values for --min");
299 : }
300 : grid_check++;
301 : }
302 9 : std::vector<std::string> gmax(cvs.size() );
303 18 : if(parseVector("--max",gmax)) {
304 4 : if(gmax.size()!=cvs.size() && gmax.size()!=0) {
305 0 : plumed_merror("not enough values for --max");
306 : }
307 4 : grid_check++;
308 : }
309 9 : std::vector<std::string> gbin(cvs.size());
310 : bool grid_has_bin;
311 : grid_has_bin=false;
312 18 : if(parseVector("--bin",gbin)) {
313 5 : if(gbin.size()!=cvs.size() && gbin.size()!=0) {
314 0 : plumed_merror("not enough values for --bin");
315 : }
316 : grid_has_bin=true;
317 : }
318 9 : std::vector<std::string> gspacing(cvs.size());
319 : bool grid_has_spacing;
320 : grid_has_spacing=false;
321 18 : if(parseVector("--spacing",gspacing)) {
322 1 : if(gspacing.size()!=cvs.size() && gspacing.size()!=0) {
323 0 : plumed_merror("not enough values for --spacing");
324 : }
325 : grid_has_spacing=true;
326 : }
327 : // allowed: no grids only bin
328 : // not allowed: partial grid definition
329 9 : plumed_massert( gmin.size()==gmax.size() && (gmin.size()==0 || gmin.size()==cvs.size() ),"you should specify --min and --max together with same number of components");
330 :
331 :
332 :
333 9 : PlumedMain plumed;
334 : std::string ss;
335 9 : unsigned nn=1;
336 : ss="setNatoms";
337 9 : plumed.cmd(ss,&nn);
338 9 : if(Communicator::initialized()) {
339 0 : plumed.cmd("setMPIComm",&pc.Get_comm());
340 : }
341 9 : plumed.cmd("init",&nn);
342 9 : std::vector <bool> isdone(cvs.size(),false);
343 26 : for(unsigned i=0; i<cvs.size(); i++) {
344 17 : if(!isdone[i]) {
345 : isdone[i]=true;
346 : std::vector<std::string> actioninput;
347 : std::vector <unsigned> inds;
348 17 : actioninput.push_back("FAKE");
349 34 : actioninput.push_back("ATOMS=1");
350 34 : actioninput.push_back("LABEL="+cvs[i][0]);
351 : std::vector<std::string> comps, periods;
352 17 : if(cvs[i].size()>1) {
353 1 : comps.push_back(cvs[i][1]);
354 1 : inds.push_back(i);
355 : }
356 17 : periods.push_back(pmin[i]);
357 17 : periods.push_back(pmax[i]);
358 25 : for(unsigned j=i+1; j<cvs.size(); j++) {
359 8 : if(cvs[i][0]==cvs[j][0] && !isdone[j]) {
360 0 : if(cvs[i].size()==1 || cvs[j].size()==1 ) {
361 0 : plumed_merror("you cannot have twice the same label and no components ");
362 : }
363 0 : if(cvs[j].size()>1) {
364 0 : comps.push_back(cvs[j][1]);
365 0 : periods.push_back(pmin[j]);
366 0 : periods.push_back(pmax[j]);
367 : isdone[j]=true;
368 0 : inds.push_back(j);
369 : }
370 : }
371 :
372 : }
373 : // drain all the components
374 : std::string addme;
375 17 : if(comps.size()>0) {
376 : addme="COMPONENTS=";
377 1 : for(unsigned ii=0; ii<comps.size()-1; ii++) {
378 0 : addme+=comps[ii]+",";
379 : }
380 : addme+=comps.back();
381 1 : actioninput.push_back(addme);
382 : }
383 : // periodicity (always explicit here)
384 : addme="PERIODIC=";
385 34 : for(unsigned j=0; j<periods.size()-1; j++) {
386 34 : addme+=periods[j]+",";
387 : }
388 : addme+=periods.back();
389 17 : actioninput.push_back(addme);
390 18 : for(unsigned j=0; j<inds.size(); j++) {
391 : unsigned jj;
392 1 : jj=inds[j];
393 1 : if(grid_check==2) {
394 : double gm;
395 : double pm;
396 1 : if(pmin[jj]!="none") {
397 1 : Tools::convert(gmin[jj],gm);
398 1 : Tools::convert(pmin[jj],pm);
399 1 : if( gm<pm ) {
400 0 : plumed_merror("Periodicity issue : GRID_MIN value ( "+gmin[jj]+" ) is less than periodicity in HILLS file in "+cvs[jj][0]+ " ( "+pmin[jj]+" ) ");
401 : }
402 : }
403 1 : if(pmax[jj]!="none") {
404 1 : Tools::convert(gmax[jj],gm);
405 1 : Tools::convert(pmax[jj],pm);
406 1 : if( gm>pm ) {
407 0 : plumed_merror("Periodicity issue : GRID_MAX value ( "+gmax[jj]+" ) is more than periodicity in HILLS file in "+cvs[jj][0]+ " ( "+pmax[jj]+" ) ");
408 : }
409 : }
410 : }
411 : }
412 :
413 : // for(unsigned i=0;i< actioninput.size();i++){
414 : // cerr<<"AA "<<actioninput[i]<<endl;
415 : // }
416 17 : plumed.readInputWords(actioninput,false);
417 34 : }
418 :
419 : }
420 9 : unsigned ncv=cvs.size();
421 : std::vector<std::string> actioninput;
422 : std::vector<std::string> idw;
423 : // check if the variables to be used are correct
424 18 : if(parseVector("--idw",idw)) {
425 4 : for(unsigned i=0; i<idw.size(); i++) {
426 : bool found=false;
427 6 : for(unsigned j=0; j<cvs.size(); j++) {
428 4 : if(cvs[j].size()>1) {
429 0 : if(idw[i]==cvs[j][0]+"."+cvs[j][1]) {
430 : found=true;
431 : }
432 : } else {
433 4 : if(idw[i]==cvs[j][0]) {
434 : found=true;
435 : }
436 : }
437 : }
438 2 : if(!found) {
439 0 : plumed_merror("variable "+idw[i]+" is not found in the bunch of cvs: revise your --idw option" );
440 : }
441 : }
442 2 : plumed_massert( idw.size()<=cvs.size(),"the number of variables to be integrated should be at most equal to the total number of cvs ");
443 : // in this case you need a beta factor!
444 : }
445 :
446 : std::string kt;
447 9 : kt=std::string("1.");// assign an arbitrary value just in case that idw.size()==cvs.size()
448 9 : if ( dohisto || idw.size()!=0 ) {
449 6 : plumed_massert(parse("--kt",kt),"if you make a dimensionality reduction (--idw) or a histogram (--histo) then you need to define --kt ");
450 : }
451 :
452 : std::string addme;
453 :
454 9 : actioninput.push_back("FUNCSUMHILLS");
455 18 : actioninput.push_back("ISCLTOOL");
456 :
457 : // set names
458 : std::string outfile;
459 18 : if(parse("--outfile",outfile)) {
460 0 : actioninput.push_back("OUTHILLS="+outfile);
461 : }
462 : std::string outhisto;
463 18 : if(parse("--outhisto",outhisto)) {
464 2 : actioninput.push_back("OUTHISTO="+outhisto);
465 : }
466 :
467 :
468 : addme="ARG=";
469 17 : for(unsigned i=0; i<(ncv-1); i++) {
470 8 : if(cvs[i].size()==1) {
471 16 : addme+=std::string(cvs[i][0])+",";
472 : } else {
473 0 : addme+=std::string(cvs[i][0])+"."+std::string(cvs[i][1])+",";
474 : }
475 : }
476 9 : if(cvs[ncv-1].size()==1) {
477 16 : addme+=std::string(cvs[ncv-1][0]);
478 : } else {
479 2 : addme+=std::string(cvs[ncv-1][0])+"."+std::string(cvs[ncv-1][1]);
480 : }
481 9 : actioninput.push_back(addme);
482 : //for(unsigned i=0;i< actioninput.size();i++){
483 : // cerr<<"AA "<<actioninput[i]<<endl;
484 : //}
485 9 : if(dohills) {
486 : addme="HILLSFILES=";
487 9 : for(unsigned i=0; i<hillsFiles.size()-1; i++) {
488 2 : addme+=hillsFiles[i]+",";
489 : }
490 : addme+=hillsFiles[hillsFiles.size()-1];
491 8 : actioninput.push_back(addme);
492 : // set the grid
493 : }
494 9 : if(grid_check==2) {
495 : addme="GRID_MAX=";
496 7 : for(unsigned i=0; i<(ncv-1); i++) {
497 6 : addme+=gmax[i]+",";
498 : }
499 : addme+=gmax[ncv-1];
500 4 : actioninput.push_back(addme);
501 : addme="GRID_MIN=";
502 7 : for(unsigned i=0; i<(ncv-1); i++) {
503 6 : addme+=gmin[i]+",";
504 : }
505 : addme+=gmin[ncv-1];
506 4 : actioninput.push_back(addme);
507 : }
508 9 : if(grid_has_bin) {
509 : addme="GRID_BIN=";
510 10 : for(unsigned i=0; i<(ncv-1); i++) {
511 10 : addme+=gbin[i]+",";
512 : }
513 : addme+=gbin[ncv-1];
514 5 : actioninput.push_back(addme);
515 : }
516 9 : if(grid_has_spacing) {
517 : addme="GRID_SPACING=";
518 1 : for(unsigned i=0; i<(ncv-1); i++) {
519 0 : addme+=gspacing[i]+",";
520 : }
521 : addme+=gspacing[ncv-1];
522 1 : actioninput.push_back(addme);
523 : }
524 : std::string stride;
525 : stride="";
526 18 : if(parse("--stride",stride)) {
527 1 : actioninput.push_back("INITSTRIDE="+stride);
528 : bool nohistory;
529 1 : parseFlag("--nohistory",nohistory);
530 1 : if(nohistory) {
531 0 : actioninput.push_back("NOHISTORY");
532 : }
533 : }
534 : bool mintozero;
535 9 : parseFlag("--mintozero",mintozero);
536 9 : if(mintozero) {
537 0 : actioninput.push_back("MINTOZERO");
538 : }
539 9 : if(idw.size()!=0) {
540 : addme="PROJ=";
541 2 : for(unsigned i=0; i<idw.size()-1; i++) {
542 0 : addme+=idw[i]+",";
543 : }
544 : addme+=idw.back();
545 2 : actioninput.push_back(addme);
546 : }
547 :
548 9 : if(dohisto) {
549 1 : if(idw.size()==0) {
550 1 : if(sigma.size()!=hcvs.size()) {
551 0 : plumed_merror("you should define as many --sigma vector as the number of collective variable used for the histogram ");
552 : }
553 : } else {
554 0 : if(idw.size()!=sigma.size()) {
555 0 : plumed_merror("you should define as many --sigma vector as the number of collective variable used for the histogram ");
556 : }
557 : }
558 : }
559 :
560 9 : if(idw.size()!=0 || dohisto) {
561 6 : actioninput.push_back("KT="+kt);
562 : }
563 9 : if(dohisto) {
564 : addme="HISTOFILES=";
565 1 : for(unsigned i=0; i<histoFiles.size()-1; i++) {
566 0 : addme+=histoFiles[i]+",";
567 : }
568 : addme+=histoFiles[histoFiles.size()-1];
569 1 : actioninput.push_back(addme);
570 :
571 : addme="HISTOSIGMA=";
572 2 : for(unsigned i=0; i<sigma.size()-1; i++) {
573 2 : addme+=sigma[i]+",";
574 : }
575 : addme+=sigma.back();
576 1 : actioninput.push_back(addme);
577 : }
578 :
579 : bool negbias;
580 9 : parseFlag("--negbias",negbias);
581 9 : if(negbias) {
582 2 : actioninput.push_back("NEGBIAS");
583 : }
584 :
585 9 : if(lowI_!=uppI_) {
586 : addme="INTERVAL=";
587 0 : addme+=lowI_+",";
588 : addme+=uppI_;
589 0 : actioninput.push_back(addme);
590 : }
591 :
592 : std::string fmt;
593 : fmt="";
594 18 : parse("--fmt",fmt);
595 9 : if(fmt!="") {
596 18 : actioninput.push_back("FMT="+fmt);
597 : }
598 :
599 :
600 : // for(unsigned i=0;i< actioninput.size();i++){
601 : // cerr<<"AA "<<actioninput[i]<<endl;
602 : // }
603 9 : plumed.readInputWords(actioninput,false);
604 : // if not a grid, then set it up automatically
605 9 : return 0;
606 27 : }
607 :
608 9 : bool CLToolSumHills::findCvsAndPeriodic(const std::string & filename,
609 : std::vector< std::vector<std::string> > &cvs,
610 : std::vector<std::string> &pmin,
611 : std::vector<std::string> &pmax,
612 : bool &multivariate,
613 : std::string &lowI_,
614 : std::string &uppI_) {
615 9 : IFile ifile;
616 9 : ifile.allowIgnoredFields();
617 : std::vector<std::string> fields;
618 9 : if(ifile.FileExist(filename)) {
619 : cvs.clear();
620 : pmin.clear();
621 : pmax.clear();
622 9 : ifile.open(filename);
623 9 : ifile.scanFieldList(fields);
624 : bool before_sigma=true;
625 121 : for(unsigned i=0; i<fields.size(); i++) {
626 : size_t pos = 0;
627 : size_t founds,foundm,foundp;
628 : //found=(fields[i].find("sigma_", pos) || fields[i].find("min_", pos) || fields[i].find("max_", pos) ) ;
629 112 : founds=fields[i].find("sigma_", pos) ;
630 112 : foundm=fields[i].find("min_", pos) ;
631 112 : foundp=fields[i].find("max_", pos) ;
632 112 : if (founds!=std::string::npos || foundm!=std::string::npos || foundp!=std::string::npos ) {
633 : before_sigma=false;
634 : }
635 : // cvs are after time and before sigmas
636 : size_t found;
637 112 : found=fields[i].find("time", pos);
638 112 : if( found==std::string::npos && before_sigma) {
639 : // separate the components
640 : size_t dot=fields[i].find_first_of('.');
641 : std::vector<std::string> ss;
642 : // this loop does not take into account repetitions
643 17 : if(dot!=std::string::npos) {
644 1 : std::string a=fields[i].substr(0,dot);
645 1 : std::string name=fields[i].substr(dot+1);
646 1 : ss.push_back(a);
647 1 : ss.push_back(name);
648 1 : cvs.push_back(ss);
649 : } else {
650 32 : cvs.emplace_back(std::vector<std::string> {fields[i]});
651 : }
652 : //std::cerr<<"found variable number "<<cvs.size()<<" : "<<cvs.back()[0]<<std::endl;
653 : //if((cvs.back()).size()!=1){
654 : // std::cerr<<"component "<<(cvs.back()).back()<<std::endl;
655 : //}
656 : // get periodicity
657 17 : pmin.push_back("none");
658 34 : pmax.push_back("none");
659 : std::string mm;
660 17 : if((cvs.back()).size()>1) {
661 2 : mm=cvs.back()[0]+"."+cvs.back()[1];
662 : } else {
663 : mm=cvs.back()[0];
664 : }
665 34 : if(ifile.FieldExist("min_"+mm)) {
666 : std::string val;
667 28 : ifile.scanField("min_"+mm,val);
668 : pmin[pmin.size()-1]=val;
669 : // std::cerr<<"found min : "<<pmin.back()<<std::endl;
670 : }
671 : //std::cerr<<"found min : "<<pmin.back()<<std::endl;
672 17 : if(ifile.FieldExist("max_"+mm)) {
673 : std::string val;
674 28 : ifile.scanField("max_"+mm,val);
675 : pmax[pmax.size()-1]=val;
676 : // std::cerr<<"found max : "<<pmax.back()<<std::endl;
677 : }
678 : //std::cerr<<"found max : "<<pmax.back()<<std::endl;
679 17 : }
680 : }
681 : // is multivariate ???
682 : std::string sss;
683 9 : multivariate=false;
684 18 : if(ifile.FieldExist("multivariate")) {
685 : ;
686 18 : ifile.scanField("multivariate",sss);
687 9 : if(sss=="true") {
688 6 : multivariate=true;
689 3 : } else if(sss=="false") {
690 3 : multivariate=false;
691 : }
692 : }
693 : // do interval?
694 18 : if(ifile.FieldExist("lower_int")) {
695 0 : ifile.scanField("lower_int",lowI_);
696 0 : ifile.scanField("upper_int",uppI_);
697 : } else {
698 : lowI_="-1.";
699 : uppI_="-1.";
700 : }
701 9 : ifile.scanField();
702 : return true;
703 : } else {
704 : return false;
705 : }
706 9 : }
707 :
708 :
709 17063 : PLUMED_REGISTER_CLTOOL(CLToolSumHills,"sum_hills")
710 :
711 :
712 :
713 : }
714 : }
|
