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 [driver](driver.md).
47 :
48 :
49 : ## 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 :
55 : ```plumed
56 : t1_1: TD_UNIFORM MINIMA=-4.0 MAXIMA=+4.0
57 : t1_2: TD_GAUSSIAN CENTER1=-2.0 SIGMA1=0.5
58 : t1: TD_PRODUCT_DISTRIBUTION DISTRIBUTIONS=t1_1,t1_2
59 :
60 : VES_OUTPUT_TARGET_DISTRIBUTION ...
61 : GRID_MIN=-4.0,-4.0
62 : GRID_MAX=+4.0,+4.0
63 : GRID_BINS=100,100
64 : TARGET_DISTRIBUTION=t1
65 : TARGETDIST_FILE=targetdist-1.data
66 : LOG_TARGETDIST_FILE=targetdist-1.log.data
67 : FMT_GRIDS=%11.6f
68 : ... VES_OUTPUT_TARGET_DISTRIBUTION
69 : ```
70 :
71 : This input should be run through the driver by using a command similar to the
72 : following one where the trajectory/configuration file configuration.gro is needed to
73 : trick the code to exit correctly.
74 :
75 : ```plumed
76 : plumed driver --plumed plumed.dat --igro configuration.gro
77 : ```
78 :
79 : */
80 : //+ENDPLUMEDOC
81 :
82 :
83 : class OutputTargetDistribution :
84 : public Action {
85 : public:
86 : explicit OutputTargetDistribution(const ActionOptions&);
87 0 : void calculate() override {}
88 0 : void apply() override {}
89 : static void registerKeywords(Keywords& keys);
90 : };
91 :
92 :
93 : PLUMED_REGISTER_ACTION(OutputTargetDistribution,"VES_OUTPUT_TARGET_DISTRIBUTION")
94 :
95 122 : void OutputTargetDistribution::registerKeywords(Keywords& keys) {
96 122 : Action::registerKeywords(keys);
97 122 : keys.add("compulsory","GRID_MIN","the lower bounds for the grid");
98 122 : keys.add("compulsory","GRID_MAX","the upper bounds for the grid");
99 122 : keys.add("compulsory","GRID_BINS","the number of bins used for the grid.");
100 122 : 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.");
101 122 : keys.add("compulsory","TARGETDIST_FILE","filename of the file for writing the target distribution");
102 122 : keys.add("optional","LOG_TARGETDIST_FILE","filename of the file for writing the log of the target distribution");
103 122 : keys.add("compulsory","TARGET_DISTRIBUTION","the target distribution to be used.");
104 122 : keys.add("optional","FMT_GRIDS","the numerical format of the target distribution grids written to file. By default it is %14.9f");
105 122 : keys.addFlag("DO_1D_PROJECTIONS",false,"Also output the one-dimensional marginal distributions for multi-dimensional target distribution.");
106 122 : }
107 :
108 120 : OutputTargetDistribution::OutputTargetDistribution(const ActionOptions&ao):
109 120 : Action(ao) {
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 120 : if(grid_periodicity.size()==0) {
131 224 : grid_periodicity.assign(nargs,"NO");
132 : }
133 :
134 120 : std::string fmt_grids="%14.9f";
135 120 : parse("FMT_GRIDS",fmt_grids);
136 :
137 120 : bool do_1d_proj = false;
138 120 : parseFlag("DO_1D_PROJECTIONS",do_1d_proj);
139 120 : if(do_1d_proj && nargs==1) {
140 0 : plumed_merror("doesn't make sense to use the DO_1D_PROJECTIONS keyword for a one-dimensional distribution");
141 : }
142 :
143 120 : plumed_massert(grid_min.size()==nargs,"mismatch between number of values given for grid parameters");
144 120 : plumed_massert(grid_max.size()==nargs,"mismatch between number of values given for grid parameters");
145 120 : plumed_massert(grid_periodicity.size()==nargs,"mismatch between number of values given for grid parameters");
146 :
147 : std::string targetdist_label;
148 120 : parse("TARGET_DISTRIBUTION",targetdist_label);
149 120 : checkRead();
150 : //
151 120 : std::vector<std::unique_ptr<Value>> arguments(nargs);
152 287 : for(unsigned int i=0; i < nargs; i++) {
153 : std::string is;
154 167 : Tools::convert(i+1,is);
155 167 : if(nargs==1) {
156 : is="";
157 : }
158 167 : arguments[i]= Tools::make_unique<Value>(nullptr,"arg"+is,false);
159 167 : if(grid_periodicity[i]=="YES") {
160 11 : arguments[i]->setDomain(grid_min[i],grid_max[i]);
161 156 : } else if(grid_periodicity[i]=="NO") {
162 156 : arguments[i]->setNotPeriodic();
163 : } else {
164 0 : plumed_merror("wrong value given in GRID_PERIODICITY, either specify YES or NO");
165 : }
166 : }
167 :
168 120 : std::string error_msg = "";
169 120 : TargetDistribution* targetdist_pntr = VesTools::getPointerFromLabel<TargetDistribution*>(targetdist_label,plumed.getActionSet(),error_msg);
170 120 : if(error_msg.size()>0) {
171 0 : plumed_merror("Error in keyword TARGET_DISTRIBUTION of "+getName()+": "+error_msg);
172 : }
173 : //
174 120 : if(targetdist_pntr->isDynamic()) {
175 0 : plumed_merror(getName() + " only works for static target distributions");
176 : }
177 120 : targetdist_pntr->setupGrids(Tools::unique2raw(arguments),grid_min,grid_max,grid_bins);
178 120 : targetdist_pntr->updateTargetDist();
179 : Grid* targetdist_grid_pntr = targetdist_pntr->getTargetDistGridPntr();
180 : Grid* log_targetdist_grid_pntr = targetdist_pntr->getLogTargetDistGridPntr();
181 :
182 :
183 120 : double sum_grid = TargetDistribution::integrateGrid(targetdist_grid_pntr);
184 120 : log.printf(" target distribution integrated over the grid: %16.12f\n",sum_grid);
185 120 : log.printf(" (%30.16e)\n",sum_grid);
186 : //
187 120 : OFile ofile;
188 120 : ofile.link(*this);
189 120 : ofile.enforceBackup();
190 120 : ofile.open(targetdist_fname);
191 : targetdist_grid_pntr->setOutputFmt(fmt_grids);
192 120 : targetdist_grid_pntr->writeToFile(ofile);
193 120 : ofile.close();
194 120 : if(log_targetdist_fname.size()>0) {
195 23 : OFile ofile2;
196 23 : ofile2.link(*this);
197 23 : ofile2.enforceBackup();
198 23 : ofile2.open(log_targetdist_fname);
199 : log_targetdist_grid_pntr->setOutputFmt(fmt_grids);
200 23 : log_targetdist_grid_pntr->writeToFile(ofile2);
201 23 : ofile2.close();
202 23 : }
203 :
204 120 : if(do_1d_proj) {
205 12 : for(unsigned int i=0; i<nargs; i++) {
206 8 : std::vector<std::string> arg1d(1);
207 8 : arg1d[0] = arguments[i]->getName();
208 8 : Grid marginal_grid = targetdist_pntr->getMarginal(arg1d);
209 : //
210 : std::string suffix;
211 8 : Tools::convert(i+1,suffix);
212 8 : suffix = "proj-" + suffix;
213 8 : std::string marginal_fname = FileBase::appendSuffix(targetdist_fname,"."+suffix);
214 : //
215 8 : OFile ofile3;
216 8 : ofile3.link(*this);
217 8 : ofile3.enforceBackup();
218 8 : ofile3.open(marginal_fname);
219 : marginal_grid.setOutputFmt(fmt_grids);
220 8 : marginal_grid.writeToFile(ofile3);
221 16 : }
222 : }
223 :
224 480 : }
225 :
226 :
227 :
228 :
229 :
230 : }
231 : }
|