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 "core/ActionPilot.h"
23 : #include "core/ActionWithValue.h"
24 : #include "core/ActionRegister.h"
25 : #include "core/PlumedMain.h"
26 : #include "core/ActionSet.h"
27 : #include "core/Atoms.h"
28 : #include "tools/IFile.h"
29 : #include <memory>
30 :
31 : namespace PLMD {
32 : namespace generic {
33 :
34 : //+PLUMEDOC GENERIC READ
35 : /*
36 : Read quantities from a colvar file.
37 :
38 : This Action can be used with driver to read in a colvar file that was generated during
39 : an MD simulation
40 :
41 : \par Description of components
42 :
43 : The READ command will read those fields that are labelled with the text string given to the
44 : VALUE keyword. It will also read in any fields that are labeled with the text string
45 : given to the VALUE keyword followed by a dot and a further string. If a single Value is read in
46 : this value can be referenced using the label of the Action. Alternatively, if multiple quantities
47 : are read in, they can be referenced elsewhere in the input by using the label for the Action
48 : followed by a dot and the character string that appeared after the dot in the title of the field.
49 :
50 : \par Examples
51 :
52 : This input reads in data from a file called input_colvar.data that was generated
53 : in a calculation that involved PLUMED. The first command reads in the data from the
54 : column headed phi1 while the second reads in the data from the column headed phi2.
55 :
56 : \plumedfile
57 : rphi1: READ FILE=input_colvar.data VALUES=phi1
58 : rphi2: READ FILE=input_colvar.data VALUES=phi2
59 : PRINT ARG=rphi1,rphi2 STRIDE=500 FILE=output_colvar.data
60 : \endplumedfile
61 :
62 : The file input_colvar.data is just a normal colvar file as shown below
63 :
64 : \auxfile{input_colvar.data}
65 : #! FIELDS time phi psi metad.bias metad.rbias metad.rct
66 : #! SET min_phi -pi
67 : #! SET max_phi pi
68 : #! SET min_psi -pi
69 : #! SET max_psi pi
70 : 0.000000 -1.2379 0.8942 0.0000 0.0000 0.0000
71 : 1.000000 -1.4839 1.0482 0.0000 0.0000 0.0089
72 : 2.000000 -1.3243 0.6055 0.0753 0.0664 0.0184
73 : \endauxfile
74 :
75 : */
76 : //+ENDPLUMEDOC
77 :
78 296 : class Read :
79 : public ActionPilot,
80 : public ActionWithValue
81 : {
82 : private:
83 : bool ignore_time;
84 : bool ignore_forces;
85 : bool cloned_file;
86 : unsigned nlinesPerStep;
87 : std::string filename;
88 : /// Unique pointer with the same scope as ifile.
89 : std::unique_ptr<IFile> ifile_ptr;
90 : /// Pointer to input file.
91 : /// It is either pointing to the content of ifile_ptr
92 : /// or to the file it is cloned from.
93 : IFile* ifile;
94 : std::vector<std::unique_ptr<Value>> readvals;
95 : public:
96 : static void registerKeywords( Keywords& keys );
97 : explicit Read(const ActionOptions&);
98 : void prepare() override;
99 90735 : void apply() override {}
100 : void calculate() override;
101 : void update() override;
102 : std::string getFilename() const;
103 : IFile* getFile();
104 : unsigned getNumberOfDerivatives() override;
105 : void turnOnDerivatives() override;
106 : };
107 :
108 7980 : PLUMED_REGISTER_ACTION(Read,"READ")
109 :
110 75 : void Read::registerKeywords(Keywords& keys) {
111 75 : Action::registerKeywords(keys);
112 75 : ActionPilot::registerKeywords(keys);
113 75 : ActionWithValue::registerKeywords(keys);
114 375 : keys.add("compulsory","STRIDE","1","the frequency with which the file should be read.");
115 375 : keys.add("compulsory","EVERY","1","only read every \\f$n\\f$th line of the colvar file. This should be used if the colvar was written more frequently than the trajectory.");
116 300 : keys.add("compulsory","VALUES","the values to read from the file");
117 300 : keys.add("compulsory","FILE","the name of the file from which to read these quantities");
118 : keys.addFlag("IGNORE_TIME",false,"ignore the time in the colvar file. When this flag is not present read will be quite strict "
119 225 : "about the start time of the simulation and the stride between frames");
120 : keys.addFlag("IGNORE_FORCES",false,"use this flag if the forces added by any bias can be safely ignored. As an example forces can be "
121 225 : "safely ignored if you are doing post processing that does not involve outputting forces");
122 150 : keys.remove("NUMERICAL_DERIVATIVES");
123 150 : keys.use("UPDATE_FROM");
124 150 : keys.use("UPDATE_UNTIL");
125 75 : ActionWithValue::useCustomisableComponents(keys);
126 75 : }
127 :
128 74 : Read::Read(const ActionOptions&ao):
129 : Action(ao),
130 : ActionPilot(ao),
131 : ActionWithValue(ao),
132 : ignore_time(false),
133 : ignore_forces(false),
134 148 : nlinesPerStep(1)
135 : {
136 : // Read the file name from the input line
137 148 : parse("FILE",filename);
138 : // Check if time is to be ignored
139 148 : parseFlag("IGNORE_TIME",ignore_time);
140 : // Check if forces are to be ignored
141 148 : parseFlag("IGNORE_FORCES",ignore_forces);
142 : // Open the file if it is not already opened
143 74 : cloned_file=false;
144 148 : std::vector<Read*> other_reads=plumed.getActionSet().select<Read*>();
145 126 : for(unsigned i=0; i<other_reads.size(); ++i) {
146 52 : if( other_reads[i]->getFilename()==filename ) {
147 26 : ifile=other_reads[i]->getFile();
148 26 : cloned_file=true;
149 : }
150 : }
151 74 : if( !cloned_file ) {
152 55 : ifile_ptr.reset(new IFile());
153 55 : ifile=ifile_ptr.get();
154 55 : if( !ifile->FileExist(filename) ) error("could not find file named " + filename);
155 55 : ifile->link(*this);
156 55 : ifile->open(filename);
157 55 : ifile->allowIgnoredFields();
158 : }
159 148 : parse("EVERY",nlinesPerStep);
160 74 : if(nlinesPerStep>1) log.printf(" only reading every %uth line of file %s\n",nlinesPerStep,filename.c_str() );
161 74 : else log.printf(" reading data from file %s\n",filename.c_str() );
162 : // Find out what we are reading
163 222 : std::vector<std::string> valread; parseVector("VALUES",valread);
164 :
165 : std::size_t dot=valread[0].find_first_of('.');
166 74 : if( valread[0].find(".")!=std::string::npos ) {
167 6 : std::string label=valread[0].substr(0,dot);
168 12 : std::string name=valread[0].substr(dot+1);
169 6 : if( name=="*" ) {
170 1 : if( valread.size()>1 ) error("all values must be from the same Action when using READ");
171 : std::vector<std::string> fieldnames;
172 1 : ifile->scanFieldList( fieldnames );
173 15 : for(unsigned i=0; i<fieldnames.size(); ++i) {
174 14 : if( fieldnames[i].substr(0,dot)==label ) {
175 9 : readvals.emplace_back(new Value(this, fieldnames[i], false) ); addComponentWithDerivatives( fieldnames[i].substr(dot+1) );
176 6 : if( ifile->FieldExist("min_" + fieldnames[i]) ) componentIsPeriodic( fieldnames[i].substr(dot+1), "-pi","pi" );
177 6 : else componentIsNotPeriodic( fieldnames[i].substr(dot+1) );
178 : }
179 1 : }
180 : } else {
181 5 : readvals.emplace_back(new Value(this, valread[0], false) ); addComponentWithDerivatives( name );
182 10 : if( ifile->FieldExist("min_" + valread[0]) ) componentIsPeriodic( valread[0].substr(dot+1), "-pi", "pi" );
183 10 : else componentIsNotPeriodic( valread[0].substr(dot+1) );
184 7 : for(unsigned i=1; i<valread.size(); ++i) {
185 2 : if( valread[i].substr(0,dot)!=label ) error("all values must be from the same Action when using READ");;
186 3 : readvals.emplace_back(new Value(this, valread[i], false) ); addComponentWithDerivatives( valread[i].substr(dot+1) );
187 2 : if( ifile->FieldExist("min_" + valread[i]) ) componentIsPeriodic( valread[i].substr(dot+1), "-pi", "pi" );
188 2 : else componentIsNotPeriodic( valread[i].substr(dot+1) );
189 : }
190 : }
191 : } else {
192 68 : if( valread.size()!=1 ) error("all values must be from the same Action when using READ");
193 68 : readvals.emplace_back(new Value(this, valread[0], false) ); addValueWithDerivatives();
194 154 : if( ifile->FieldExist("min_" + valread[0]) ) setPeriodic( "-pi", "pi" );
195 59 : else setNotPeriodic();
196 68 : log.printf(" reading value %s and storing as %s\n",valread[0].c_str(),getLabel().c_str() );
197 : }
198 74 : checkRead();
199 74 : }
200 :
201 26 : std::string Read::getFilename() const {
202 26 : return filename;
203 : }
204 :
205 26 : IFile* Read::getFile() {
206 26 : return ifile;
207 : }
208 :
209 0 : unsigned Read::getNumberOfDerivatives() {
210 0 : return 0;
211 : }
212 :
213 20 : void Read::turnOnDerivatives() {
214 20 : if( !ignore_forces ) error("cannot calculate derivatives for colvars that are read in from a file. If you are postprocessing and "
215 0 : "these forces do not matter add the flag IGNORE_FORCES to all READ actions");
216 20 : }
217 :
218 90735 : void Read::prepare() {
219 90735 : if( !cloned_file ) {
220 : double du_time;
221 88890 : if( !ifile->scanField("time",du_time) ) {
222 0 : error("Reached end of file " + filename + " before end of trajectory");
223 88890 : } else if( fabs( du_time-getTime() )>plumed.getAtoms().getTimeStep() && !ignore_time ) {
224 0 : std::string str_dutime,str_ptime; Tools::convert(du_time,str_dutime); Tools::convert(getTime(),str_ptime);
225 0 : error("mismatched times in colvar files : colvar time=" + str_dutime + " plumed time=" + str_ptime + ". Add IGNORE_TIME to ignore error.");
226 : }
227 : }
228 90735 : }
229 :
230 90735 : void Read::calculate() {
231 : std::string smin, smax;
232 365326 : for(unsigned i=0; i<readvals.size(); ++i) {
233 : // .get returns the raw pointer
234 : // ->get calls the Value::get() method
235 91928 : ifile->scanField( readvals[i].get() );
236 91928 : getPntrToComponent(i)->set( readvals[i]->get() );
237 91928 : if( readvals[i]->isPeriodic() ) {
238 3309 : readvals[i]->getDomain( smin, smax );
239 3309 : getPntrToComponent(i)->setDomain( smin, smax );
240 : }
241 : }
242 90735 : }
243 :
244 90735 : void Read::update() {
245 90735 : if( !cloned_file ) {
246 44445 : for(unsigned i=0; i<nlinesPerStep; ++i) {
247 44445 : ifile->scanField(); double du_time;
248 132226 : if( plumed.getAtoms().getNatoms()==0 && !ifile->scanField("time",du_time) ) plumed.stop();
249 : }
250 : }
251 90735 : }
252 :
253 : }
254 5874 : }
|
