Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2017 of Haochuan Chen (excluding colvar_UIestimator.h)
3 : Copyright (c) 2017 of Haohao Fu (colvar_UIestimator.h)
4 :
5 : This program is free software: you can redistribute it and/or modify
6 : it under the terms of the GNU Lesser General Public License as published
7 : by the Free Software Foundation, either version 3 of the License, or
8 : (at your option) any later version.
9 :
10 : This program is distributed in the hope that it will be useful,
11 : but WITHOUT ANY WARRANTY; without even the implied warranty of
12 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 : GNU Lesser General Public License for more details.
14 :
15 : You should have received a copy of the GNU Lesser General Public License
16 : along with this program. If not, see <http://www.gnu.org/licenses/>.
17 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
18 : #ifdef __PLUMED_HAS_BOOST_SERIALIZATION
19 : #include "cltools/CLTool.h"
20 : #include "cltools/CLToolRegister.h"
21 : #include "config/Config.h"
22 : #include "core/ActionRegister.h"
23 : #include "DRR.h"
24 : #include "tools/Tools.h"
25 : #include "tools/Units.h"
26 : #include <boost/archive/binary_iarchive.hpp>
27 : #include <boost/archive/binary_oarchive.hpp>
28 : #include <boost/serialization/vector.hpp>
29 : #include <cstdlib>
30 : #include <fstream>
31 : #include <iostream>
32 : #include <string>
33 :
34 : using namespace PLMD;
35 : using namespace cltools;
36 :
37 : namespace PLMD {
38 : namespace drr {
39 :
40 : //+PLUMEDOC EABFMOD_TOOLS drr_tool
41 : /*
42 : - Extract .grad and .count files from the binary output .drrstate
43 : - Merge windows
44 :
45 : \par Examples
46 :
47 : The following command will extract .grad and .count files.
48 : \verbatim
49 : plumed drr_tool --extract eabf.drrstate
50 : \endverbatim
51 :
52 : The following command will merge windows of two .drrstate file, and output the
53 : .grad and .count files.
54 : \verbatim
55 : plumed drr_tool --merge win1.drrstate,win2.drrstate
56 : \endverbatim
57 :
58 : After getting the .grad and .count file, you can do numerical integration by
59 : using abf_integrate tool from
60 : https://github.com/Colvars/colvars/tree/master/colvartools
61 : \verbatim
62 : abf_integrate eabf.czar.grad
63 : \endverbatim
64 : \note
65 : The abf_integrate in colvartools is in kcal/mol, so it may be better to use --units kcal/mol when running drr_tool
66 :
67 : */
68 : //+ENDPLUMEDOC
69 :
70 9 : class drrtool : public CLTool {
71 : public:
72 : static void registerKeywords(Keywords &keys);
73 : explicit drrtool(const CLToolOptions &co);
74 : int main(FILE *in, FILE *out, Communicator &pc);
75 : void extractdrr(const std::vector<std::string> &filename);
76 : void mergewindows(const std::vector<std::string> &filename);
77 0 : std::string description() const { return "Extract or merge the drrstate files."; }
78 :
79 : private:
80 : bool verbosity;
81 : Units units;
82 : const std::string suffix{".drrstate"};
83 : };
84 :
85 6455 : PLUMED_REGISTER_CLTOOL(drrtool, "drr_tool")
86 :
87 1613 : void drrtool::registerKeywords(Keywords &keys) {
88 1613 : CLTool::registerKeywords(keys);
89 6452 : keys.add("optional", "--extract", "Extract drrstate file(s)");
90 6452 : keys.add("optional", "--merge", "Merge eABF windows");
91 8065 : keys.add("compulsory","--units","kj/mol","the units of energy can be kj/mol, kcal/mol, j/mol, eV or the conversion factor from kj/mol");
92 4839 : keys.addFlag("-v", false, "Verbose output");
93 1613 : }
94 :
95 3 : drrtool::drrtool(const CLToolOptions &co) : CLTool(co) {
96 3 : inputdata = commandline;
97 3 : verbosity = false;
98 3 : }
99 :
100 3 : int drrtool::main(FILE *in, FILE *out, Communicator &pc) {
101 6 : parseFlag("-v", verbosity);
102 3 : std::vector<std::string> stateFilesToExtract;
103 : std::string unitname;
104 6 : parse("--units",unitname);
105 3 : units.setEnergy( unitname );
106 6 : bool doextract = parseVector("--extract", stateFilesToExtract);
107 3 : if (doextract) {
108 2 : extractdrr(stateFilesToExtract);
109 : }
110 3 : std::vector<std::string> stateFilesToMerge;
111 6 : bool domerge = parseVector("--merge", stateFilesToMerge);
112 3 : if (domerge) {
113 1 : mergewindows(stateFilesToMerge);
114 : }
115 3 : return 0;
116 : }
117 :
118 2 : void drrtool::extractdrr(const std::vector<std::string> &filename) {
119 6 : #pragma omp parallel for
120 4 : for (size_t j = 0; j < filename.size(); ++j) {
121 4 : std::ifstream in;
122 4 : in.open(filename[j]);
123 : boost::archive::binary_iarchive ia(in);
124 : long long int step;
125 : std::vector<double> fict;
126 : std::vector<double> vfict;
127 : std::vector<double> vfict_laststep;
128 : std::vector<double> ffict;
129 : ABF abfgrid;
130 : CZAR czarestimator;
131 : ia >> step >> fict >> vfict >> vfict_laststep >> ffict >> abfgrid >>
132 : czarestimator;
133 2 : in.close();
134 2 : abfgrid.setOutputUnit(units.getEnergy());
135 : czarestimator.setOutputUnit(units.getEnergy());
136 2 : if (verbosity) {
137 1 : std::cout << "Output units factor: " << units.getEnergy() << '\n';
138 : std::cout << "Dumping information of extended variables..." << '\n';
139 1 : std::cout << "Step: " << step << '\n';
140 3 : for (size_t i = 0; i < fict.size(); ++i) {
141 1 : std::cout << "Dimension[" << i + 1 << "]:\n"
142 1 : << " Coordinate: " << fict[i] << '\n'
143 1 : << " Velocity: " << vfict[i] << '\n'
144 1 : << " Velocity(laststep): " << vfict_laststep[i] << '\n'
145 1 : << " Force: " << ffict[i] << '\n';
146 : }
147 : std::cout << "Dumping counts and gradients from grids..." << '\n';
148 : }
149 2 : std::string outputname(filename[j]);
150 4 : outputname = outputname.substr(0, outputname.length() - suffix.length());
151 2 : if (verbosity)
152 : std::cout << "Writing ABF(naive) estimator files..." << '\n';
153 2 : abfgrid.writeAll(outputname);
154 2 : if (verbosity)
155 : std::cout << "Writing CZAR estimator files..." << '\n';
156 2 : czarestimator.writeAll(outputname);
157 : }
158 2 : }
159 :
160 1 : void drrtool::mergewindows(const std::vector<std::string> &filename) {
161 1 : if (filename.size() < 2) {
162 : std::cerr << "ERROR! You need at least two .drrstate file to merge windows!" << std::endl;
163 0 : std::abort();
164 : }
165 : // Read grid into abfs and czars;
166 1 : std::vector<ABF> abfs;
167 1 : std::vector<CZAR> czars;
168 3 : for (auto it_fn = filename.begin(); it_fn != filename.end(); ++it_fn) {
169 4 : std::ifstream in;
170 2 : in.open((*it_fn));
171 : boost::archive::binary_iarchive ia(in);
172 : long long int step;
173 : std::vector<double> fict;
174 : std::vector<double> vfict;
175 : std::vector<double> vfict_laststep;
176 : std::vector<double> ffict;
177 : ABF abfgrid;
178 : CZAR czarestimator;
179 : ia >> step >> fict >> vfict >> vfict_laststep >> ffict >> abfgrid >>
180 : czarestimator;
181 2 : abfgrid.setOutputUnit(units.getEnergy());
182 : czarestimator.setOutputUnit(units.getEnergy());
183 2 : abfs.push_back(abfgrid);
184 2 : czars.push_back(czarestimator);
185 2 : in.close();
186 : }
187 1 : CZAR cmerged = CZAR::mergewindow(czars[0], czars[1]);
188 1 : ABF amerged = ABF::mergewindow(abfs[0], abfs[1]);
189 1 : for (size_t i = 2; i < czars.size(); ++i) {
190 0 : cmerged = CZAR::mergewindow(cmerged, czars[i]);
191 0 : amerged = ABF::mergewindow(amerged, abfs[i]);
192 : }
193 : // Generate new file name for merged grad and count
194 2 : std::vector<std::string> tmp_name = filename;
195 6 : std::transform(std::begin(tmp_name), std::end(tmp_name), std::begin(tmp_name), [&](std::string s) {return s.substr(0, s.find(suffix));});
196 7 : std::string mergename = std::accumulate(std::begin(tmp_name), std::end(tmp_name), std::string(""), [](std::string a, std::string b) {return a + b + "+";});
197 2 : mergename = mergename.substr(0, mergename.size() - 1);
198 1 : cmerged.writeAll(mergename);
199 1 : amerged.writeAll(mergename);
200 1 : }
201 :
202 : } // End of namespace
203 4839 : }
204 :
205 : #endif
|
