Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2012-2019 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 "PathMSDBase.h"
23 : #include "core/PlumedMain.h"
24 :
25 : using namespace std;
26 :
27 : namespace PLMD {
28 : namespace colvar {
29 :
30 : //+PLUMEDOC COLVAR PROPERTYMAP
31 : /*
32 : Calculate generic property maps.
33 :
34 : This Colvar calculates the property maps according to the work of Spiwok \cite Spiwok:2011ce.
35 :
36 :
37 : Basically it calculates
38 : \f{eqnarray}
39 : X=\frac{\sum_i X_i*\exp(-\lambda D_i(x))}{\sum_i \exp(-\lambda D_i(x))} \\
40 : Y=\frac{\sum_i Y_i*\exp(-\lambda D_i(x))}{\sum_i \exp(-\lambda D_i(x))} \\
41 : \cdots\\
42 : zzz=-\frac{1}{\lambda}\log(\sum_i \exp(-\lambda D_i(x)))
43 : \f}
44 :
45 : where the parameters \f$X_i\f$ and \f$Y_i\f$ are provided in the input pdb (allv.pdb in this case) and
46 : \f$D_i(x)\f$ is the MSD after optimal alignment calculated on the pdb frames you input (see Kearsley).
47 :
48 : \warning
49 : The molecule used for \ref PROPERTYMAP calculation should be whole (both atoms used in alignment and in displacement calculation).
50 : In case it is broken by the host MD code, please use \ref WHOLEMOLECULES to reconstruct it before \ref PROPERTYMAP calculation.
51 :
52 : \par Examples
53 :
54 : \plumedfile
55 : p3: PROPERTYMAP REFERENCE=../../trajectories/path_msd/allv.pdb PROPERTY=X,Y LAMBDA=69087 NEIGH_SIZE=8 NEIGH_STRIDE=4
56 : PRINT ARG=p3.X,p3.Y,p3.zzz STRIDE=1 FILE=colvar FMT=%8.4f
57 : \endplumedfile
58 :
59 : note that NEIGH_STRIDE=4 NEIGH_SIZE=8 control the neighborlist parameter (optional but
60 : recommended for performance) and states that the neighbor list will be calculated every 4
61 : timesteps and consider only the closest 8 member to the actual md snapshots.
62 :
63 : In this case the input line instructs plumed to look for two properties X and Y with attached values in the REMARK
64 : line of the reference pdb (Note: No spaces from X and = and 1 !!!!).
65 : e.g.
66 :
67 : \verbatim
68 : REMARK X=1 Y=2
69 : ATOM 1 CL ALA 1 -3.171 0.295 2.045 1.00 1.00
70 : ATOM 5 CLP ALA 1 -1.819 -0.143 1.679 1.00 1.00
71 : .......
72 : END
73 : REMARK X=2 Y=3
74 : ATOM 1 CL ALA 1 -3.175 0.365 2.024 1.00 1.00
75 : ATOM 5 CLP ALA 1 -1.814 -0.106 1.685 1.00 1.00
76 : ....
77 : END
78 : \endverbatim
79 :
80 : \note
81 : The implementation of this collective variable and of \ref PATHMSD
82 : is shared, as well as most input options.
83 :
84 : */
85 : //+ENDPLUMEDOC
86 :
87 22 : class PropertyMap : public PathMSDBase {
88 : public:
89 : explicit PropertyMap(const ActionOptions&);
90 : static void registerKeywords(Keywords& keys);
91 : };
92 :
93 6463 : PLUMED_REGISTER_ACTION(PropertyMap,"PROPERTYMAP")
94 :
95 12 : void PropertyMap::registerKeywords(Keywords& keys) {
96 12 : PathMSDBase::registerKeywords(keys);
97 48 : keys.add("compulsory","PROPERTY","the property to be used in the indexing: this goes in the REMARK field of the reference");
98 12 : ActionWithValue::useCustomisableComponents(keys);
99 48 : keys.addOutputComponent("zzz","default","the minimum distance from the reference points");
100 12 : }
101 :
102 11 : PropertyMap::PropertyMap(const ActionOptions&ao):
103 : Action(ao),
104 11 : PathMSDBase(ao)
105 : {
106 : // this is the only additional keyword needed
107 22 : parseVector("PROPERTY",labels);
108 11 : checkRead();
109 11 : log<<" Bibliography "
110 33 : <<plumed.cite("Spiwok V, Kralova B J. Chem. Phys. 135, 224504 (2011)")
111 22 : <<"\n";
112 11 : if(labels.size()==0) {
113 : char buf[500];
114 : sprintf(buf,"Need to specify PROPERTY with this action\n");
115 0 : plumed_merror(buf);
116 : } else {
117 88 : for(unsigned i=0; i<labels.size(); i++) {
118 44 : log<<" found custom propety to be found in the REMARK line: "<<labels[i].c_str()<<"\n";
119 66 : addComponentWithDerivatives(labels[i]); componentIsNotPeriodic(labels[i]);
120 : }
121 : // add distance anyhow
122 33 : addComponentWithDerivatives("zzz"); componentIsNotPeriodic("zzz");
123 : //reparse the REMARK field and pick the index
124 1408 : for(unsigned i=0; i<pdbv.size(); i++) {
125 924 : vector<std::string> myv(pdbv[i].getRemark());
126 : // now look for X=1.34555 Y=5.6677
127 : vector<double> labelvals;
128 3696 : for(unsigned j=0; j<labels.size(); j++) {
129 : double val;
130 924 : if(Tools::parse(myv,labels[j],val)) {labelvals.push_back(val);}
131 : else {
132 : char buf[500];
133 : sprintf(buf,"PROPERTY LABEL \" %s \" NOT FOUND IN REMARK FOR FRAME %u \n",labels[j].c_str(),i);
134 0 : plumed_merror(buf);
135 : };
136 : }
137 462 : indexvec.push_back(labelvals);
138 : }
139 : }
140 11 : requestAtoms(pdbv[0].getAtomNumbers());
141 :
142 11 : }
143 :
144 : }
145 4839 : }
146 :
147 :
|
