Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2016-2021 The VES code team
3 : (see the PEOPLE-VES file at the root of this folder for a list of names)
4 :
5 : See http://www.ves-code.org for more information.
6 :
7 : This file is part of VES code module.
8 :
9 : The VES code module 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 : The VES code module 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 the VES code module. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 :
23 : #include "TargetDistribution.h"
24 :
25 : #include "GridIntegrationWeights.h"
26 : #include "VesTools.h"
27 :
28 : #include "core/ActionRegister.h"
29 : #include "core/PlumedMain.h"
30 : #include "core/Value.h"
31 : #include "tools/File.h"
32 : #include "tools/Grid.h"
33 :
34 :
35 : namespace PLMD {
36 : namespace ves {
37 :
38 : //+PLUMEDOC VES_UTILS VES_OUTPUT_TARGET_DISTRIBUTION
39 : /*
40 : Output target distribution to file.
41 :
42 : This action can be used to output target distributions to a grid file,
43 : for example to see how they look like before using them in a VES bias.
44 : This action only support static target distributions.
45 :
46 : This action is normally used through the \ref driver.
47 :
48 :
49 : \par Examples
50 :
51 : In the following input we define a target distribution that is uniform for
52 : argument 1 and a Gaussian for argument 2 and then output it to a file
53 : called targetdist-1.data.
54 : \plumedfile
55 : t1_1: TD_UNIFORM MINIMA=-4.0 MAXIMA=+4.0
56 : t1_2: TD_GAUSSIAN CENTER1=-2.0 SIGMA1=0.5
57 : t1: TD_PRODUCT_DISTRIBUTION DISTRIBUTIONS=t1_1,t1_2
58 :
59 : VES_OUTPUT_TARGET_DISTRIBUTION ...
60 : GRID_MIN=-4.0,-4.0
61 : GRID_MAX=+4.0,+4.0
62 : GRID_BINS=100,100
63 : TARGET_DISTRIBUTION=t1
64 : TARGETDIST_FILE=targetdist-1.data
65 : LOG_TARGETDIST_FILE=targetdist-1.log.data
66 : FMT_GRIDS=%11.6f
67 : ... VES_OUTPUT_TARGET_DISTRIBUTION
68 : \endplumedfile
69 :
70 : This input should be run through the driver by using a command similar to the
71 : following one where the trajectory/configuration file configuration.gro is needed to
72 : trick the code to exit correctly.
73 : \verbatim
74 : plumed driver --plumed plumed.dat --igro configuration.gro
75 : \endverbatim
76 :
77 : */
78 : //+ENDPLUMEDOC
79 :
80 :
81 : class OutputTargetDistribution :
82 : public Action
83 : {
84 : public:
85 : explicit OutputTargetDistribution(const ActionOptions&);
86 0 : void calculate() override {}
87 0 : void apply() override {}
88 : static void registerKeywords(Keywords& keys);
89 : };
90 :
91 :
92 10539 : PLUMED_REGISTER_ACTION(OutputTargetDistribution,"VES_OUTPUT_TARGET_DISTRIBUTION")
93 :
94 121 : void OutputTargetDistribution::registerKeywords(Keywords& keys) {
95 121 : Action::registerKeywords(keys);
96 242 : keys.add("compulsory","GRID_MIN","the lower bounds for the grid");
97 242 : keys.add("compulsory","GRID_MAX","the upper bounds for the grid");
98 242 : keys.add("compulsory","GRID_BINS","the number of bins used for the grid.");
99 363 : keys.add("optional","GRID_PERIODICITY","specify if the individual arguments should be made periodic (YES) or not (NO). By default all arguments are taken as not periodic.");
100 242 : keys.add("compulsory","TARGETDIST_FILE","filename of the file for writing the target distribution");
101 363 : keys.add("optional","LOG_TARGETDIST_FILE","filename of the file for writing the log of the target distribution");
102 363 : keys.add("compulsory","TARGET_DISTRIBUTION","the target distribution to be used.");
103 242 : keys.add("optional","FMT_GRIDS","the numerical format of the target distribution grids written to file. By default it is %14.9f");
104 242 : keys.addFlag("DO_1D_PROJECTIONS",false,"Also output the one-dimensional marginal distributions for multi-dimensional target distribution.");
105 121 : }
106 :
107 120 : OutputTargetDistribution::OutputTargetDistribution(const ActionOptions&ao):
108 120 : Action(ao)
109 : {
110 :
111 : std::string targetdist_fname;
112 240 : parse("TARGETDIST_FILE",targetdist_fname);
113 : std::string log_targetdist_fname;
114 120 : parse("LOG_TARGETDIST_FILE",log_targetdist_fname);
115 120 : if(targetdist_fname==log_targetdist_fname) {
116 0 : plumed_merror("error in " + getName() + ":TARGETDIST_FILE and LOG_TARGETDIST_FILE cannot be the same");
117 : }
118 :
119 : std::vector<unsigned int> grid_bins;
120 240 : parseVector("GRID_BINS",grid_bins);
121 120 : unsigned int nargs = grid_bins.size();
122 :
123 120 : std::vector<std::string> grid_min(nargs);
124 120 : parseVector("GRID_MIN",grid_min);
125 120 : std::vector<std::string> grid_max(nargs);
126 120 : parseVector("GRID_MAX",grid_max);
127 :
128 120 : std::vector<std::string> grid_periodicity(nargs);
129 240 : parseVector("GRID_PERIODICITY",grid_periodicity);
130 232 : if(grid_periodicity.size()==0) {grid_periodicity.assign(nargs,"NO");}
131 :
132 120 : std::string fmt_grids="%14.9f";
133 120 : parse("FMT_GRIDS",fmt_grids);
134 :
135 120 : bool do_1d_proj = false;
136 120 : parseFlag("DO_1D_PROJECTIONS",do_1d_proj);
137 120 : if(do_1d_proj && nargs==1) {
138 0 : plumed_merror("doesn't make sense to use the DO_1D_PROJECTIONS keyword for a one-dimensional distribution");
139 : }
140 :
141 120 : plumed_massert(grid_min.size()==nargs,"mismatch between number of values given for grid parameters");
142 120 : plumed_massert(grid_max.size()==nargs,"mismatch between number of values given for grid parameters");
143 120 : plumed_massert(grid_periodicity.size()==nargs,"mismatch between number of values given for grid parameters");
144 :
145 : std::string targetdist_label;
146 120 : parse("TARGET_DISTRIBUTION",targetdist_label);
147 120 : checkRead();
148 : //
149 120 : std::vector<std::unique_ptr<Value>> arguments(nargs);
150 287 : for(unsigned int i=0; i < nargs; i++) {
151 167 : std::string is; Tools::convert(i+1,is);
152 167 : if(nargs==1) {is="";}
153 167 : arguments[i]= Tools::make_unique<Value>(nullptr,"arg"+is,false);
154 167 : if(grid_periodicity[i]=="YES") {
155 11 : arguments[i]->setDomain(grid_min[i],grid_max[i]);
156 : }
157 156 : else if(grid_periodicity[i]=="NO") {
158 156 : arguments[i]->setNotPeriodic();
159 : }
160 : else {
161 0 : plumed_merror("wrong value given in GRID_PERIODICITY, either specify YES or NO");
162 : }
163 : }
164 :
165 120 : std::string error_msg = "";
166 120 : TargetDistribution* targetdist_pntr = VesTools::getPointerFromLabel<TargetDistribution*>(targetdist_label,plumed.getActionSet(),error_msg);
167 120 : if(error_msg.size()>0) {plumed_merror("Error in keyword TARGET_DISTRIBUTION of "+getName()+": "+error_msg);}
168 : //
169 120 : if(targetdist_pntr->isDynamic()) {
170 0 : plumed_merror(getName() + " only works for static target distributions");
171 : }
172 120 : targetdist_pntr->setupGrids(Tools::unique2raw(arguments),grid_min,grid_max,grid_bins);
173 120 : targetdist_pntr->updateTargetDist();
174 : Grid* targetdist_grid_pntr = targetdist_pntr->getTargetDistGridPntr();
175 : Grid* log_targetdist_grid_pntr = targetdist_pntr->getLogTargetDistGridPntr();
176 :
177 :
178 120 : double sum_grid = TargetDistribution::integrateGrid(targetdist_grid_pntr);
179 120 : log.printf(" target distribution integrated over the grid: %16.12f\n",sum_grid);
180 120 : log.printf(" (%30.16e)\n",sum_grid);
181 : //
182 120 : OFile ofile;
183 120 : ofile.link(*this);
184 120 : ofile.enforceBackup();
185 120 : ofile.open(targetdist_fname);
186 : targetdist_grid_pntr->setOutputFmt(fmt_grids);
187 120 : targetdist_grid_pntr->writeToFile(ofile);
188 120 : ofile.close();
189 120 : if(log_targetdist_fname.size()>0) {
190 23 : OFile ofile2;
191 23 : ofile2.link(*this);
192 23 : ofile2.enforceBackup();
193 23 : ofile2.open(log_targetdist_fname);
194 : log_targetdist_grid_pntr->setOutputFmt(fmt_grids);
195 23 : log_targetdist_grid_pntr->writeToFile(ofile2);
196 23 : ofile2.close();
197 23 : }
198 :
199 120 : if(do_1d_proj) {
200 12 : for(unsigned int i=0; i<nargs; i++) {
201 8 : std::vector<std::string> arg1d(1);
202 8 : arg1d[0] = arguments[i]->getName();
203 8 : Grid marginal_grid = targetdist_pntr->getMarginal(arg1d);
204 : //
205 : std::string suffix;
206 8 : Tools::convert(i+1,suffix);
207 8 : suffix = "proj-" + suffix;
208 8 : std::string marginal_fname = FileBase::appendSuffix(targetdist_fname,"."+suffix);
209 : //
210 8 : OFile ofile3;
211 8 : ofile3.link(*this);
212 8 : ofile3.enforceBackup();
213 8 : ofile3.open(marginal_fname);
214 : marginal_grid.setOutputFmt(fmt_grids);
215 8 : marginal_grid.writeToFile(ofile3);
216 16 : }
217 : }
218 :
219 480 : }
220 :
221 :
222 :
223 :
224 :
225 : }
226 : }
|
