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 "Path.h"
23 : #include "core/ActionRegister.h"
24 : #include "core/ActionWithValue.h"
25 : #include "core/ActionWithArguments.h"
26 : #include "core/PlumedMain.h"
27 : #include "core/ActionSet.h"
28 : #include "tools/PDB.h"
29 :
30 : //+PLUMEDOC COLVAR PATH
31 : /*
32 : Path collective variables with a more flexible framework for the distance metric being used.
33 :
34 : The Path Collective Variables developed by Branduardi and that are described in the first paper that is cited below allow one
35 : to compute the progress along a high-dimensional path and the distance from the high-dimensional
36 : path. The progress along the path (s) is computed using:
37 :
38 : $$
39 : s = \frac{ \sum_{i=1}^N i \exp( -\lambda R[X - X_i] ) }{ \sum_{i=1}^N \exp( -\lambda R[X - X_i] ) }
40 : $$
41 :
42 : while the distance from the path (z) is measured using:
43 :
44 : $$
45 : z = -\frac{1}{\lambda} \ln\left[ \sum_{i=1}^N \exp( -\lambda R[X - X_i] ) \right]
46 : $$
47 :
48 : In these expressions $N$ high-dimensional frames ($X_i$) are used to describe the path in the high-dimensional
49 : space. The two expressions above are then functions of the distances from each of the high-dimensional frames $R[X - X_i]$.
50 : Within PLUMED there are multiple ways to define the distance from a high-dimensional configuration. You could calculate
51 : the RMSD distance or you could calculate the amount by which a set of collective variables change. As such this shortcut
52 : of the path CV allows one to compute the distances from the paths in a variety of different ways. This is as opposed to
53 : the alternative implementation of path ([PATHMSD](PATHMSD.md)) which is a bit faster but which only allows one to use the [RMSD](RMSD.md) distance.
54 :
55 : ## Examples
56 :
57 : In the example below the path is defined using RMSD distance from frames.
58 :
59 : ```plumed
60 : #SETTINGS INPUTFILES=regtest/trajectories/path_msd/all.pdb
61 : p1: PATH REFERENCE=regtest/trajectories/path_msd/all.pdb TYPE=OPTIMAL LAMBDA=69087
62 : PRINT ARG=p1.spath,p1.zpath STRIDE=1 FILE=colvar
63 : ```
64 :
65 : In the example below the path is defined using the values of two torsional angles (t1 and t2). You can specify these arguments in the PLUMED input
66 : as illustrated below:
67 :
68 : ```plumed
69 : #SETTINGS INPUTFILES=regtest/mapping/rt-tpath/epath.pdb
70 : t1: TORSION ATOMS=5,7,9,15
71 : t2: TORSION ATOMS=7,9,15,17
72 : pp: PATH ARG=t1,t2 REFERENCE=regtest/mapping/rt-tpath/epath.pdb LAMBDA=1.0
73 : PRINT ARG=pp.* FILE=colvar
74 : ```
75 :
76 : However, this is not strictly necessary as PLUMED can also get the names of the arguments from the input pdb file directly as illustrated in the following
77 : input.
78 :
79 : ```plumed
80 : #SETTINGS INPUTFILES=regtest/mapping/rt-tpath/epath.pdb
81 : t1: TORSION ATOMS=5,7,9,15
82 : t2: TORSION ATOMS=7,9,15,17
83 : pp: PATH REFERENCE=regtest/mapping/rt-tpath/epath.pdb LAMBDA=1.0
84 : PRINT ARG=pp.* FILE=colvar
85 : ```
86 :
87 : If you look at the pdb file for the input above you can see that the remarks tell PLUMED the labels that are being used to define the position in the
88 : high dimensional space and the values that these arguments have at each point on the path.
89 :
90 : ## Controlling what is calculated
91 :
92 : If you only want to calculate the $s$ coordinate you can use the NOZPATH flag as shown below:
93 :
94 : ```plumed
95 : #SETTINGS INPUTFILES=regtest/trajectories/path_msd/all.pdb
96 : p1: PATH REFERENCE=regtest/trajectories/path_msd/all.pdb TYPE=OPTIMAL NOZPATH LAMBDA=69087
97 : PRINT ARG=p1.spath STRIDE=1 FILE=colvar
98 : ```
99 :
100 : Similarly, if you only want to calculate the $z$ coordinate you use the NOSPATH flag as shown below
101 :
102 : ```plumed
103 : #SETTINGS INPUTFILES=regtest/trajectories/path_msd/all.pdb
104 : p1: PATH REFERENCE=regtest/trajectories/path_msd/all.pdb TYPE=OPTIMAL NOSPATH LAMBDA=69087
105 : PRINT ARG=p1.zpath STRIDE=1 FILE=colvar
106 : ```
107 :
108 : You can also tell PLUMED to calculate the distance along and distance from the path using the formulas that are given in the documentation for the
109 : [GPATH](GPATH.md) shortcut as shown below.
110 :
111 : ```plumed
112 : #SETTINGS INPUTFILES=regtest/trajectories/path_msd/all.pdb
113 : p1: PATH REFERENCE=regtest/trajectories/path_msd/all.pdb TYPE=OPTIMAL GPATH LAMBDA=69087
114 : PRINT ARG=p1.spath,p1.zpath,p1.gspath,p1.gzpath STRIDE=1 FILE=colvar
115 : ```
116 :
117 : The PATH command in this input calculates four quantities:
118 :
119 : - p1.spath - the s coordinate evaluated using the expression above
120 : - p1.zpath - the z coordinate evluated using the expression above
121 : - p1.gspath - the s coordinate evaluated using the expression in the documentation for [GPATH](GPATH.md)
122 : - p1.gzpath - the z coordinate evaluated using the expression in the documentation for [GPATH](GPATH.md)
123 :
124 : ## Paths using normalized Euclidean distances
125 :
126 : If you are defining the positions of the frames in the path using arguments you may wish to calculate the distances between the instantaneous values of the
127 : arguments and the values of these arguments on the path using the [NORMALIZED_EUCLIDEAN_DISTANCE](NORMALIZED_EUCLIDEAN_DISTANCE.md) rather than the
128 : [EUCLIDEAN_DISTANCE](EUCLIDEAN_DISTANCE.md). In other words, instead of calculating the distance between the vectors coordinates of coordinates $x$ and $y$ using:
129 :
130 : $$
131 : D^2 = \sum_{i=1}^N (x_i - y_i)^2
132 : $$
133 :
134 : where the sum runs over the number of arguments, you can calculate this distance as:
135 :
136 : $$
137 : D^2 = \sum_{i=1}^N C_i (x_i - y_i)^2
138 : $$
139 :
140 : where the $C_i$ are a set of coefficients that are specified using the COEFFICIENTS keyword as shown below.
141 :
142 : ```plumed
143 : #SETTINGS INPUTFILES=regtest/mapping/rt-tpath/epath.pdb
144 : t1: TORSION ATOMS=5,7,9,15
145 : t2: TORSION ATOMS=7,9,15,17
146 : pp: PATH ...
147 : ARG=t1,t2 COEFFICIENTS=0.2,0.3
148 : REFERENCE=regtest/mapping/rt-tpath/epath.pdb
149 : LAMBDA=1.0
150 : ...
151 : PRINT ARG=pp.* FILE=colvar
152 : ```
153 :
154 : */
155 : //+ENDPLUMEDOC
156 :
157 : //+PLUMEDOC COLVAR GPROPERTYMAP
158 : /*
159 : Property maps but with a more flexible framework for the distance metric being used.
160 :
161 : This colvar calculates a property map using the formalism developed by Spiwok that is referenced in the second paper cited below.
162 : In essence if you have the value of some property, $X_i$, that it takes at a set of high-dimensional
163 : positions then you calculate the value of the property at some arbitrary point in the high-dimensional space
164 : using:
165 :
166 : $$
167 : X=\frac{\sum_i X_i*\exp(-\lambda D_i(x))}{\sum_i \exp(-\lambda D_i(x))}
168 : $$
169 :
170 : In these expressions the value of the property $X$ is given at the position of $N$ high-dimensional frames ($X_i$). The distances, $D_i$, between
171 : these frames and our instantaneous position are then used to compute the value of the property at the instantaneous position in the high dimensional space.
172 : As illustrated in the examples below, this implementation of the property map allows you to calculate these distances in various different ways.
173 : This is as opposed to the alternative implementation [PROPERTYMAP](PROPERTYMAP.md) which is a bit faster but which only allows one to use the RMSD distance.
174 :
175 : ## Examples
176 :
177 : The input shown below can be used to calculate the interpolated values of two properties called X and Y based on the values
178 : that these properties take at a set of reference configurations and using the formula above. For this input the distances
179 : between the reference configurations and the instantaneous configurations are calculated using the OPTIMAL metric that is
180 : discussed at length in the manual pages on [RMSD](RMSD.md).
181 :
182 : ```plumed
183 : #SETTINGS INPUTFILES=regtest/trajectories/path_msd/allv.pdb
184 : p2: GPROPERTYMAP REFERENCE=regtest/trajectories/path_msd/allv.pdb PROPERTY=X,Y LAMBDA=69087
185 : PRINT ARG=p2_X,p2_Y,p2_zpath STRIDE=1 FILE=colvar
186 : ```
187 :
188 : Notice that the REFERENCE file here gives the values of the properties at each of the points of interest.
189 : In this second input the value of the property at each point of interest is also given in the REFERENCE file.
190 : Here as well the REFERENCE file also tells us that each of the reference points is defined based on the values of the
191 : two torsions `t1` and `t2`.
192 :
193 : ```plumed
194 : #SETTINGS INPUTFILES=regtest/mapping/rt-tpath/epath.pdb
195 : t1: TORSION ATOMS=5,7,9,15
196 : t2: TORSION ATOMS=7,9,15,17
197 : pp: GPROPERTYMAP ...
198 : ARG=t1,t2 REFERENCE=regtest/mapping/rt-tpath/epath.pdb
199 : PROPERTY=X LAMBDA=1.0
200 : ...
201 : PRINT ARG=pp.* FILE=colvar
202 : ```
203 :
204 : Notice that you can use the ARG keyword to specify the arguments are used in the definition of the PATH. However, you do not have to
205 : use this keyword. The following input also works as PLUMED gets the input arguments from the file that contains the definition of the
206 : path.
207 :
208 : ```plumed
209 : #SETTINGS INPUTFILES=regtest/mapping/rt-tpath/epath.pdb
210 : t1: TORSION ATOMS=5,7,9,15
211 : t2: TORSION ATOMS=7,9,15,17
212 : pp: GPROPERTYMAP ...
213 : REFERENCE=regtest/mapping/rt-tpath/epath.pdb
214 : PROPERTY=X LAMBDA=1.0
215 : ...
216 : PRINT ARG=pp.* FILE=colvar
217 :
218 : ## Controlling what is calculated
219 :
220 : If you dont want to bother calculating the z coordinate you can use the NOZPATH flag as shown below:
221 :
222 : ```plumed
223 : #SETTINGS INPUTFILES=regtest/trajectories/path_msd/allv.pdb
224 : p2: GPROPERTYMAP ...
225 : REFERENCE=regtest/trajectories/path_msd/allv.pdb
226 : PROPERTY=X,Y LAMBDA=69087 NOZPATH
227 : ...
228 : PRINT ARG=p2_X,p2_Y STRIDE=1 FILE=colvar
229 : ```
230 :
231 : ## Property maps using normalized Euclidean distances
232 :
233 : If you are defining the positions of the frames in the path using arguments you may wish to calculate the distances between the instantaneous values of the
234 : arguments and the values of these arguments on the path using the [NORMALIZED_EUCLIDEAN_DISTANCE](NORMALIZED_EUCLIDEAN_DISTANCE.md) rather than the
235 : [EUCLIDEAN_DISTANCE](EUCLIDEAN_DISTANCE.md). In other words, instead of calculating the distance between the vectors coordinates of coordinates $x$ and $y$ using:
236 :
237 : $$
238 : D^2 = \sum_{i=1}^N (x_i - y_i)^2
239 : $$
240 :
241 : where the sum runs over the number of arguments, you can calculate this distance as:
242 :
243 : $$
244 : D^2 = \sum_{i=1}^N C_i (x_i - y_i)^2
245 : $$
246 :
247 : where the $C_i$ are a set of coefficients that are specified using the COEFFICIENTS keyword as shown below.
248 :
249 : ```plumed
250 : #SETTINGS INPUTFILES=regtest/mapping/rt-tpath/epath.pdb
251 : t1: TORSION ATOMS=5,7,9,15
252 : t2: TORSION ATOMS=7,9,15,17
253 : pp: GPROPERTYMAP ...
254 : REFERENCE=regtest/mapping/rt-tpath/epath.pdb
255 : PROPERTY=X COEFFICIENTS=0.2,0.3 LAMBDA=1.0
256 : ...
257 : PRINT ARG=pp.* FILE=colvar
258 : ```
259 :
260 : */
261 : //+ENDPLUMEDOC
262 :
263 : namespace PLMD {
264 : namespace mapping {
265 :
266 : PLUMED_REGISTER_ACTION(Path,"PATH")
267 : PLUMED_REGISTER_ACTION(Path,"GPROPERTYMAP")
268 :
269 26 : void Path::registerKeywords( Keywords& keys ) {
270 26 : ActionShortcut::registerKeywords( keys );
271 26 : Path::registerInputFileKeywords( keys );
272 26 : if( keys.getDisplayName()=="GPROPERTYMAP" ) {
273 7 : keys.add("optional","PROPERTY","the property to be used in the index. This should be in the REMARK of the reference");
274 : }
275 26 : keys.add("optional","LAMBDA","the lambda parameter is needed for smoothing, is in the units of plumed. It is required unless the spath and zpath are not calculated (NOSPATH and NOZPATH)");
276 26 : keys.addFlag("NOZPATH",false,"do not calculate the zpath CV");
277 26 : keys.add("optional","COEFFICIENTS","the coefficients of the displacements along each argument that should be used when calculating the euclidean distance");
278 26 : if( keys.getDisplayName()=="PATH" ) {
279 19 : keys.addFlag("NOSPATH",false,"do not calculate the spath CV");
280 19 : keys.addFlag("GPATH",false,"calculate the geometric path");
281 38 : keys.addOutputComponent("gspath","GPATH","scalar","the position along the path calculated using the geometric formula");
282 38 : keys.addOutputComponent("gzpath","GPATH","scalar","the distance from the path calculated using the geometric formula");
283 : }
284 52 : keys.addOutputComponent("spath","default","scalar","the position along the path calculated");
285 52 : keys.addOutputComponent("zpath","default","scalar","the distance from the path calculated");
286 26 : keys.addDOI("10.1063/1.2432340");
287 26 : keys.addDOI("10.1063/1.3660208");
288 26 : }
289 :
290 59 : void Path::registerInputFileKeywords( Keywords& keys ) {
291 59 : keys.add("compulsory","REFERENCE","a pdb file containing the set of reference configurations");
292 59 : keys.add("compulsory","TYPE","OPTIMAL-FAST","the manner in which distances are calculated. More information on the different "
293 : "metrics that are available in PLUMED can be found in the documentation for RMSD");
294 118 : keys.addInputKeyword("optional","ARG","scalar","the list of arguments you would like to use in your definition of the path");
295 59 : keys.needsAction("DRMSD");
296 59 : keys.needsAction("RMSD");
297 59 : keys.needsAction("LOWEST");
298 59 : keys.needsAction("GPATH");
299 59 : keys.needsAction("EUCLIDEAN_DISTANCE");
300 59 : keys.needsAction("CUSTOM");
301 59 : keys.needsAction("SUM");
302 59 : keys.needsAction("COMBINE");
303 59 : keys.needsAction("NORMALIZED_EUCLIDEAN_DISTANCE");
304 59 : keys.needsAction("PDB2CONSTANT");
305 59 : keys.needsAction("CONSTANT");
306 59 : }
307 :
308 14 : Path::Path( const ActionOptions& ao ):
309 : Action(ao),
310 14 : ActionShortcut(ao) {
311 : // Now create all other parts of the calculation
312 : std::string reference;
313 14 : parse("REFERENCE",reference);
314 14 : bool nospath=false, nozpath=false, gpath=false;
315 14 : if( getName()=="PATH" ) {
316 11 : parseFlag("NOSPATH",nospath);
317 22 : parseFlag("GPATH",gpath);
318 : }
319 30 : parseFlag("NOZPATH",nozpath);
320 : // LAMBDA is only used by spath and zpath: it is not needed with e.g. GPATH NOSPATH NOZPATH
321 : std::string lambda;
322 28 : parse("LAMBDA",lambda);
323 14 : if( lambda.length()==0 && !(nospath && nozpath) ) {
324 2 : error("LAMBDA is required to calculate the spath and zpath (use NOSPATH and NOZPATH if you do not need them)");
325 : }
326 12 : if( gpath ) {
327 2 : readInputLine( getShortcutLabel() + "_gpath: GPATH REFERENCE=" + reference );
328 4 : readInputLine( getShortcutLabel() + "_gspath: COMBINE ARG=" + getShortcutLabel() + "_gpath.s PERIODIC=NO");
329 8 : readInputLine( getShortcutLabel() + "_gzpath: COMBINE ARG=" + getShortcutLabel() + "_gpath.z PERIODIC=NO");
330 : }
331 : // the reference frames are read with these keywords (even if only the geometric path is calculated,
332 : // so that they are not reported as unknown keywords)
333 : std::string type, reference_data;
334 12 : std::vector<std::string> argnames;
335 12 : parseVector("ARG",argnames);
336 12 : parse("TYPE",type);
337 12 : if( nospath && nozpath ) {
338 : return;
339 : }
340 : // Setup the properties
341 11 : std::vector<std::string> properties, pnames;
342 11 : if( getName()=="PATH") {
343 8 : properties.resize(1);
344 : } else {
345 6 : parseVector("PROPERTY",pnames);
346 3 : properties.resize( pnames.size() );
347 : }
348 11 : if( type.find("DRMSD")!=std::string::npos ) {
349 3 : readInputLine( getShortcutLabel() + "_data: DRMSD SQUARED TYPE=" + type + " REFERENCE=" + reference );
350 : } else {
351 10 : readInputFrames( reference, type, argnames, false, this, reference_data );
352 : }
353 : // Find the shortest distance to the frames
354 22 : readInputLine( getShortcutLabel() + "_mindist: LOWEST ARG=" + getShortcutLabel() + "_data");
355 : // Now create MATHEVAL object to compute exponential functions
356 44 : readInputLine( getShortcutLabel() + "_weights: CUSTOM ARG=" + getShortcutLabel() + "_data," + getShortcutLabel() + "_mindist FUNC=exp(-(x-y)*" + lambda + ") PERIODIC=NO" );
357 : // Create denominator
358 22 : readInputLine( getShortcutLabel() + "_denom: SUM ARG=" + getShortcutLabel() + "_weights PERIODIC=NO");
359 : // Now compte zpath variable
360 11 : if( !nozpath ) {
361 44 : readInputLine( getShortcutLabel() + "_z: CUSTOM ARG=" + getShortcutLabel() + "_denom," + getShortcutLabel() + "_mindist FUNC=y-log(x)/" + lambda + " PERIODIC=NO");
362 33 : readInputLine( getShortcutLabel() + "_zpath: COMBINE ARG=" + getShortcutLabel() + "_z PERIODIC=NO");
363 : }
364 : // Now get coefficients for properies for spath
365 11 : readPropertyInformation( pnames, getShortcutLabel(), reference, this );
366 : // Now create COMBINE objects to compute numerator of path
367 25 : for(unsigned i=0; i<properties.size(); ++i) {
368 14 : if( pnames.size()>0 ) {
369 24 : readInputLine( pnames[i] + "_numer_prod: CUSTOM ARG=" + getShortcutLabel() + "_weights," + pnames[i] + "_ref FUNC=x*y PERIODIC=NO");
370 12 : readInputLine( pnames[i] + "_numer: SUM ARG=" + pnames[i] + "_numer_prod PERIODIC=NO");
371 30 : readInputLine( getShortcutLabel() + "_" + pnames[i] + ": CUSTOM ARG=" + pnames[i] + "_numer," + getShortcutLabel() + "_denom FUNC=x/y PERIODIC=NO");
372 8 : } else if( !nospath ) {
373 24 : readInputLine( getShortcutLabel() + "_s_prod: CUSTOM ARG=" + getShortcutLabel() + "_weights," + getShortcutLabel() + "_ind FUNC=x*y PERIODIC=NO");
374 16 : readInputLine( getShortcutLabel() + "_numer: SUM ARG=" + getShortcutLabel() + "_s_prod PERIODIC=NO");
375 24 : readInputLine( getShortcutLabel() + "_s: CUSTOM ARG=" + getShortcutLabel() + "_numer," + getShortcutLabel() + "_denom FUNC=x/y PERIODIC=NO");
376 24 : readInputLine( getShortcutLabel() + "_spath: COMBINE ARG=" + getShortcutLabel() + "_s PERIODIC=NO");
377 : }
378 : }
379 16 : }
380 :
381 19 : std::string Path::fixArgumentName( const std::string& argin ) {
382 19 : std::string argout = argin;
383 19 : std::size_t dot=argin.find(".");
384 19 : if( dot!=std::string::npos ) {
385 12 : argout = argin.substr(0,dot) + "_" + argin.substr(dot+1);
386 : }
387 19 : return argout;
388 : }
389 :
390 20 : void Path::readInputFrames( const std::string& reference, const std::string& type, std::vector<std::string>& argnames, const bool& displacements, ActionShortcut* action, std::string& reference_data ) {
391 20 : FILE* fp=std::fopen(reference.c_str(),"r");
392 20 : PDB pdb;
393 20 : if(!fp) {
394 0 : action->error("could not open reference file " + reference );
395 : }
396 20 : bool do_read=pdb.readFromFilepointer(fp,false,0.1);
397 20 : if( !do_read ) {
398 0 : action->error("missing file " + reference );
399 : }
400 20 : if( pdb.getPositions().size()!=0 && argnames.size()==0 ) {
401 15 : reference_data = action->getShortcutLabel() + "_data_ref";
402 15 : if( displacements ) {
403 21 : action->readInputLine( action->getShortcutLabel() + "_data: RMSD DISPLACEMENT SQUARED REFERENCE=" + reference + " TYPE=" + type );
404 : } else {
405 24 : action->readInputLine( action->getShortcutLabel() + "_data: RMSD SQUARED REFERENCE=" + reference + " TYPE=" + type );
406 : }
407 5 : } else if( pdb.getPositions().size()!=0 ) {
408 0 : reference_data = action->getShortcutLabel() + "_atomdata_ref";
409 0 : if( displacements ) {
410 0 : action->readInputLine( action->getShortcutLabel() + "_atomdata: RMSD DISPLACEMENT SQUARED REFERENCE=" + reference + " TYPE=" + type );
411 : } else {
412 0 : action->readInputLine( action->getShortcutLabel() + "_atomdata: RMSD SQUARED REFERENCE=" + reference + " TYPE=" + type );
413 : }
414 5 : } else if( argnames.size()==0 ) {
415 0 : argnames.resize( pdb.getArgumentNames().size() );
416 0 : for(unsigned i=0; i<argnames.size(); ++i) {
417 0 : argnames[i] = pdb.getArgumentNames()[i];
418 : }
419 : }
420 20 : std::vector<Value*> theargs;
421 20 : if( argnames.size()>0 ) {
422 5 : ActionWithArguments::interpretArgumentList( argnames, action->plumed.getActionSet(), action, theargs );
423 : }
424 :
425 20 : if( theargs.size()>0 ) {
426 : std::string instargs, refargs;
427 20 : for(unsigned i=0; i<theargs.size(); ++i) {
428 15 : std::string iargn = fixArgumentName( theargs[i]->getName() );
429 45 : action->readInputLine( action->getShortcutLabel() + "_ref_" + iargn + ": PDB2CONSTANT REFERENCE=" + reference + " ARG=" + theargs[i]->getName() );
430 15 : if( i==0 ) {
431 5 : instargs=" ARG2=" + theargs[i]->getName();
432 15 : refargs=" ARG1=" + action->getShortcutLabel() + "_ref_" + iargn;
433 : } else {
434 10 : instargs +="," + theargs[i]->getName();
435 30 : refargs +="," + action->getShortcutLabel() + "_ref_" + iargn;
436 : }
437 15 : if( pdb.getPositions().size()==0 && i==0 ) {
438 10 : reference_data = action->getShortcutLabel() + "_ref_" + iargn;
439 : } else {
440 30 : reference_data += "," + action->getShortcutLabel() + "_ref_" + iargn;
441 : }
442 : }
443 5 : std::string comname="EUCLIDEAN_DISTANCE SQUARED";
444 : std::string coeffstr;
445 5 : if( action->keywords.exists("COEFFICIENTS") ) {
446 4 : action->parse("COEFFICIENTS",coeffstr);
447 : }
448 5 : if( coeffstr.length()>0 ) {
449 1 : if( displacements ) {
450 0 : action->error("cannot use COEFFICIENTS arguments with GEOMETRIC PATH");
451 : }
452 1 : action->readInputLine( action->getShortcutLabel() + "_coeff: CONSTANT VALUES=" + coeffstr );
453 2 : action->readInputLine( action->getShortcutLabel() + "_coeff2: CUSTOM ARG=" + action->getShortcutLabel() + "_coeff FUNC=x*x PERIODIC=NO");
454 3 : comname = "NORMALIZED_EUCLIDEAN_DISTANCE SQUARED METRIC=" + action->getShortcutLabel() + "_coeff2";
455 4 : } else if( displacements ) {
456 : comname = "DISPLACEMENT";
457 : }
458 :
459 5 : if( pdb.getPositions().size()==0 ) {
460 5 : if( displacements ) {
461 2 : action->readInputLine( action->getShortcutLabel() + "_dataP: " + comname + instargs + refargs );
462 6 : action->readInputLine( action->getShortcutLabel() + "_data: CUSTOM ARG=" + action->getShortcutLabel() + "_dataP FUNC=-x PERIODIC=NO");
463 : } else {
464 6 : action->readInputLine( action->getShortcutLabel() + "_data: " + comname + instargs + refargs );
465 : }
466 : } else {
467 0 : action->readInputLine( action->getShortcutLabel() + "_argdata: " + comname + instargs + refargs );
468 : }
469 : }
470 20 : }
471 :
472 21 : void Path::readPropertyInformation( const std::vector<std::string>& pnames, const std::string& lab, const std::string& refname, ActionShortcut* action ) {
473 21 : if( pnames.size()>0 ) {
474 15 : for(unsigned i=0; i<pnames.size(); ++i) {
475 20 : action->readInputLine( pnames[i] + ": CONSTANT VALUE=1");
476 30 : action->readInputLine( pnames[i] + "_ref: PDB2CONSTANT REFERENCE=" + refname + " ARG=" + pnames[i] );
477 : }
478 : } else {
479 16 : ActionWithValue* av=action->plumed.getActionSet().selectWithLabel<ActionWithValue*>( lab + "_data" );
480 16 : unsigned nfram = av->copyOutput(0)->getShape()[0];
481 16 : std::string indices = "VALUES=1";
482 616 : for(unsigned i=1; i<nfram; ++i) {
483 : std::string num;
484 600 : Tools::convert( i+1, num );
485 600 : indices += "," + num;
486 : }
487 32 : action->readInputLine( lab + "_ind: CONSTANT " + indices );
488 : }
489 621 : }
490 :
491 : }
492 : }
|
